2017-03-11 4 views
0

ジョン・フォックスとの共同研究により、effects CRANパッケージから機能Effectで引数xlevelsのヘルプは言う:部分的な残差のあるプロットをエフェクトする:x値が100の値で評価されるのはなぜですか?

(...) If partial residuals are 
computed, then the focal predictor that is to appear on the 
horizontal axis of an effect plot is evaluated at 100 equally 
spaced values along its full range, and, by default, other 
numeric predictors are evaluated at the quantiles specified 
in the ‘quantiles’ argument, unless their values are given 
explicitly in ‘xlevels’. 

が、私は理解していない: 1)焦点の予測はのペアでxlevelsの関連性があるもの連続変数と離散変数(私はそれが関係ないと思っていたでしょう)。 2)ほとんどの場合、コンディショニングプロットの数が約2〜4の場合、連続予測子のペアについてxlevelsの値が100の理由 3)どのようにxlevelsが横座標の焦点予測子に影響しますか(それはそうではないようです)。

1)モデルに連続的な予測子と係数のみが含まれている場合は、xlevelsは役割を果たさず、部分的な残差は予測値を実際の値で評価することによって得られます。この場合、上記のテキストは適用されないようです。次のコードは、

 library(effects) 
    library(ggplot2) 

    x11() 
    x11() 

    set.seed(123) 
    n <- 6 ## with small n it is easier to compare the plots 
    x1 <- rnorm(n, mean = 15) 
    xf <- rep(c(0, 1), c(n/2, n - n/2)) 
    dd <- data.frame(y = x1 + x1 * 4 * xf + rnorm(n, 0, sd = 0.2), 
        x1 = x1, xf = factor(xf)) 
    mi1 <- lm(y ~ x1 * xf, data = dd) 


    plot(Effect(c("x1", "xf"), mi1, partial.residuals = TRUE), nrow = 1) 

    plot(Effect(c("x1", "xf"), mi1, partial.residuals = TRUE, 
      xlevels = list(x1 = c(1, 9, 37)))) 

    ############# "By hand", step by step 
    r1 <- resid(mi1) 
    ## Individual beta_j * x_j terms to add 
    add_x1 <- x1 * coefficients(mi1)["x1"] 
    add_x1_f1 <- model.matrix(mi1)[, "x1:xf1"] * coefficients(mi1)["x1:xf1"] 

    ## Partial residuals 
    partial_residuals <- r1 + add_x1 + add_x1_f1 

    ## For convenience for ggplot2 
    pd <- data.frame(x1 = dd$x1, xf = dd$xf, 
        partial_residuals = partial_residuals, 
        fitted = fitted(mi1)) 

    dev.set(dev.next()) 
    ## Identical to the Effect plot. 
    ggplot(data = pd, aes(x = x1, y = partial_residuals)) + 
      geom_smooth(method = "lm") + geom_point() + 
      facet_wrap(~ xf, nrow = 1) 

2「手で」得られた部分残差とEffectプロットを比較し、3)。モデルに連続する2つの予測変数が含まれている場合は、xlevelsを変更することでパネルの数を変更できますが、デフォルトで100の値を評価する理由はわかりません。私が正しく理解すれば、ここで実装されていることは、John FoxとSanford Weisbergが、例えばVisualizing Lack of Fit in Complex Regression Models: Adding Partial Residuals to Effect Displaysで議論した考えです。たとえば、スライド34では、「部分的な残差を得るために、水平軸を定義する焦点予測子以外のエフェクトで予測子を丸めて、スライスされた値に丸めます。これは私が理解していますが、スライスは水平軸(コンディショニングパネルの焦点予測器のみ)の焦点予測器には影響しません。次のコードは、予測値の1つのみを手動でスライスして "手"で得られた部分残差と再度比較します。また、水平軸の焦点予測子のxlevelsを変更しても差がないように見えます。他の予測因子:

 set.seed(123) 
    n <- 12 ## with small n easier to see 
    x1 <- rnorm(n, mean = 15) 
    x2 <- rep(c(1, 2, 8, 9), length.out = n) ## Simpler if few values 
    dd2 <- data.frame(y = x1 + 2 * x1 * x2 + rnorm(n), 
         x1 = x1, x2 = x2) 
    mi12 <- lm(y ~ x1 * x2, data = dd2) 
    summary(mi12) 


    ## Residuals + Individual beta_j * x_j terms to add 
    tt <- c("x1", "x1:x2") 
    pr2 <- resid(mi12) + model.matrix(mi12)[, tt] %*% coefficients(mi12)[tt] 

    ## Partial residuals with slicing of x2: 
    ## Evaluate x2 only at 1 and 9; all values 1 and 2 are set to 1 
    ## and all 8 and 9 are set to 9. I think this is similar to 
    ## "fitted <- y[good][closest(trans(x.fit), x[good])]" 
    ## in line 317 in plot-methods.R? 

    x22 <- dd2$x2 
    x22[x22 == 2] <- 1 
    x22[x22 == 8] <- 9 

    mm2 <- cbind(x1 = dd2$x1, "x1:x2" = dd2$x1 * x22) 
    pr2b <- resid(mi12) + mm2 %*% coefficients(mi12)[tt] 

    ## For convenience for ggplot2 
    pd2 <- data.frame(x1 = dd2$x1, 
         x2 = dd2$x2, 
         x22 = x22, 
         pr2 = pr2, 
         pr2b = pr2b) 


    plot(Effect(c("x1", "x2"), mi12, partial.residuals = TRUE)) 
    dev.set(dev.next()) 
    ggplot(data = pd2, aes(x = x1, y = pr2)) + 
     geom_smooth(method = "lm") + geom_point() + 
     facet_wrap(~ factor(x2, levels = c(8, 9, 1, 2)), nrow = 2) 


    ## Only at two values of x2 
    plot(Effect(c("x1", "x2"), mi12, partial.residuals = TRUE, 
       xlevels = list(x2 = c(1, 9), x1 = 5))) 
    dev.set(dev.next()) 
    ggplot(data = pd2, aes(x = x1, y = pr2b)) + 
     geom_smooth(method = "lm") + geom_point() + 
     facet_wrap(~ factor(x22), nrow = 1) 

    ## other values for x1 make no difference 
    plot(Effect(c("x1", "x2"), mi12, partial.residuals = TRUE, 
        xlevels = list(x2 = c(1, 9), x1 = 1))) 

ファイルplot-methods.Rのライン317の周りに、コードを見て、上記の例で述べたように、(機能plot.eff)私たちはfitted <- y[good][closest(trans(x.fit), x[good])]を参照してください。 xlevelsを変更すると、例1)のxのサイズには影響がなく、x2の場合はxlevelsを変更すると、例2-3の場合にのみ違いがあります。これは、「コンディショニングパネル内のプレディクタを、離散プレディクタで連続的に不要な場合は、スライスされた値のセットに設定する」というアイデアと一貫しているようです。

答えて

0

ジョン・フォックスは答えを提供する:

を数値予測因子は、各スライスは、ディスプレイに パネルに対応するので、これは値の 比較的少数で行われる「スライス」されている場合。数値のプレディクタ (グラフの水平軸に表示される に表示されるプレディクタ)に100個の値を使用すると、比較的小さな点で適合値 を計算し、最も近い数値を に設定することができますこれは計算を高速化する。

関連する問題