2016-04-29 12 views
1

ためにロジスティック回帰からの係数を抽出、このリンクでいくつかのデータを検索します:私は巣の成功を予測するロジスティック回帰を持っているループ

https://www.dropbox.com/s/okp2iudnace6fha/data1.csv?dl=0

すべての私の説明変数は線形トレンドを得るために連続しています。 (LD)

日付を敷設

(... 0,1,2,3など)

年:私は、巣の生存の季節変動の経時変化があるかどうか分析したいですNestAge

は、これは私のモデルである:

glm(survive~LD+yr+yr^2+LD:yr+LD:yr^2+NestAge,family=binomial(link=logexp(data1$exposure)), data=data1) 

これは、私が使用していたリンクの暴露である:

library(MASS) 
    logexp <- function(exposure = 1) 
    { 
     linkfun <- function(mu) qlogis(mu^(1/exposure)) 
     ## FIXME: is there some trick we can play here to allow 
     ## evaluation in the context of the 'data' argument? 
     linkinv <- function(eta) plogis(eta)^exposure 
     mu.eta <- function(eta) exposure * plogis(eta)^(exposure-1) * 
     .Call(stats:::C_logit_mu_eta, eta, PACKAGE = "stats") 
     valideta <- function(eta) TRUE 
     link <- paste("logexp(", deparse(substitute(exposure)), ")", 
        sep="") 
     structure(list(linkfun = linkfun, linkinv = linkinv, 
        mu.eta = mu.eta, valideta = valideta, 
        name = link), 
       class = "link-glm") 
    } 

係数を取得するには、ループを作成して抽出しますが、それは機能しません。

a<-as.matrix(coef(model)) 
intercept<-a[1,] 
slope<-a[2,] 
for (i in 1:6) { 
    i<-as.numeric(i) 
    sub<-subset(data1,data1$yr==i) 
    g<- intercept + slope*sub$yr[i] 
    dsr <-exp(g)/ (1+ exp(g)) 
} 

解決してもらえますか?前もって感謝します。

+1

YRAがラインにまで言及されているもの: 'サブ<-subset(DATA1、DATA1の$のYRA == I)'?これはdata1の変数の1つではありません。 –

+1

また、データセット内の変数のどれも26レベルの要素で構成されていない場合、なぜ26回ループしていますか?ネスト、年齢、ネスト、年齢、暴露、生存、運命、LD、NestAgeのレベルを与える長さを与える 'sapply(names(data1)、function(x)table(data1 [、x])%>%length) 6,699,6,30,2,33,33。 –

+0

それは本当です!申し訳ありませんが、私のコードを編集しました。完全なデータセットであるため、26でした(私はドロップボックスにデータの例を残しておきたい)。 – MSS

答えて

0

それが動作するようになりました:

dsr=0 
for (i in 1:6) { 
    i<-as.numeric(i) 
    sub<-subset(data1,data1$yr==i) 
    g<- intercept + slope*sub$yr[i] 
    dsr [i]<-plogis(g) 
} 
関連する問題