4

多くの X 変数と Y 変数 (500 x 500 など) があります。次の小さなデータ:

yvars <- data.frame (Yv1 = rnorm(100, 5, 3), Y2 = rnorm (100, 6, 4),
  Yv3 = rnorm (100, 14, 3))
xvars <- data.frame (Xv1 = sample (c(1,0, -1), 100, replace = T),
 X2 = sample (c(1,0, -1), 100, replace = T), 
 Xv3 = sample (c(1,0, -1), 100, replace = T), 
 D = sample (c(1,0, -1), 100, replace = T))

p値を抽出して、次のような行列を作成したい:

     Yv1    Y2     Yv3   
Xv1
X2
Xv3
D

プロセスをループする私の試みは次のとおりです。

prob = NULL
   anova.pmat <- function (x) {
            mydata <- data.frame(yvar = yvars[, x], xvars)
            for (i in seq(length(xvars))) {
              prob[[i]] <- anova(lm(yvar ~ mydata[, i + 1],
              data = mydata))$`Pr(>F)`[1]
              }
              }
    sapply (yvars,anova.pmat)
    Error in .subset(x, j) : only 0's may be mixed with negative subscripts
What could be the solution ?

編集：

最初の Y 変数の場合:

最初の Y 変数の場合:

prob <- NULL
 mydata <- data.frame(yvar = yvars[, 1], xvars)
for (i in seq(length(xvars))) {
              prob[[i]] <- anova(lm(yvar ~ mydata[, i + 1],
              data = mydata))$`Pr(>F)`[1]
              }

prob 
[1] 0.4995179 0.4067040 0.4181571 0.6291167

もう一度編集します。

for (j in seq(length (yvars))){
        prob <- NULL
        mydata <- data.frame(yvar = yvars[, j], xvars)
         for (i in seq(length(xvars))) {
                  prob[[i]] <- anova(lm(yvar ~ mydata[, i + 1],
                  data = mydata))$`Pr(>F)`[1]
                  }
}

Gives the same result as above !!!
4

2 に答える 2

4

plyrこれは、データフレームの列を (リストとして扱い) をループし、適切な p 値を返し、それを行列に配置xvarsするアプローチです。yvars行/列の名前を追加するのは余分です。

library("plyr")

probs <- laply(xvars, function(x) {
    laply(yvars, function(y) {
        anova(lm(y~x))$`Pr(>F)`[1]
    })
})
rownames(probs) <- names(xvars)
colnames(probs) <- names(yvars)
于 2012-07-10T20:18:06.607 に答える
1

combnこれは、Y変数とX変数のすべての組み合わせを生成してテストし（ は使用できません）、それぞれの場合に線形モデルを実行する1つのソリューションです。

dfrm <- data.frame(y=gl(ncol(yvars), ncol(xvars), labels=names(yvars)),
                   x=gl(ncol(xvars), 1, labels=names(xvars)), pval=NA)
## little helper function to create formula on the fly
fm <- function(x) as.formula(paste(unlist(x), collapse="~"))
## merge both datasets
full.df <- cbind.data.frame(yvars, xvars)
## apply our LM row-wise
dfrm$pval <- apply(dfrm[,1:2], 1, 
                   function(x) anova(lm(fm(x), full.df))$`Pr(>F)`[1])
## arrange everything in a rectangular matrix of p-values
res <- matrix(dfrm$pval, nc=3, dimnames=list(levels(dfrm$x), levels(dfrm$y)))

補足:高次元のデータセットでは、線形回帰の p 値を計算するために QR 分解に頼るのは時間がかかります。ペアごとの比較ごとにピアソン線形相関の行列を計算し、自由度が定義されている関係 F = ν a r 2 /(1-r 2 ) を使用して r 統計量を Fisher-Snedecor F に変換する方が簡単です。 as ν a =(n-2)-#{(x i =NA),(y i =NA)} (つまり、(n-2) から対ごとの欠損値の数を引いたもの - 欠損値がない場合値の場合、この式は回帰における通常の係数 R 2です)。

于 2012-07-10T20:10:59.943 に答える