2016-06-23 6 views
1

R에있는 MCMCglmm 패키지에서 실행중인 혼합 효과 모델이 있습니다. 탐색에 관심이있는 여러 가지 예상 변수가 있습니다. 모델에서 정보가없는 예측 변수를 제외시키기 위해 모델 단순화를 수행하려면 어떻게해야합니까? 다른 상황에서MCMCglmm에 대한 모델 단순화는 어떻게 수행합니까?

, 나는 그들을 제거하기 위해 update() 다음에, 내가 관심 오전 예측을 모두 포함하는 "전체"모델로 시작하고, 가치가없는 변수를 찾기 위해 drop1()를 사용하여 즉, 회귀 단계적 뒤로 사용했다 모델.

이러한 함수 중 어느 것도 MCMCglmm 개체에서 작동하지 않습니다. 최선의 대안은 무엇입니까? drop1()이 작동하지 않는 이유 중 하나는 MCMCglmm 개체에 AIC 값이 없다는 것입니다. 대신 DIC 값을 비슷한 방식으로 사용할 수 있습니까? 나는 DIC 값에 대해 많이 알지 못하지만 아래 예제에서 AIC처럼 행동하지 않는 것 같습니다. 모델에서 변수를 제거한 후 DIC는 100을 넘었습니다.

update()의 경우, 분명히 모델을 처음부터 다시 지정할 수는 있지만 꽤 긴 함수 호출을 복사/붙여 넣기가 포함되며 상호 작용을 포함하는 더 긴 모델 수식을 사용하면 훨씬 더 복잡해집니다. 나는 더 간단한 대안이 있는지 궁금해했다.

MCMCglmm 개체에서 작동하는 함수가 작성되지 않은 이유가 있는지 궁금합니다. 어떤 식 으로든 잘못된 방법으로 시도하려고합니다.

부분적으로 발명 한 데이터를 사용하여 재현 예 :

set.seed(1234) 
library(MCMCglmm) 

data(bird.families) 
n <- Ntip(bird.families) 

# Create some dummy variables 
d <- data.frame(taxon = bird.families$tip.label, 
       X1 = rnorm(n), 
       X2 = rnorm(n), 
       X3 = sample(c("A", "B", "C"), n, replace = T), 
       X4 = sample(c("A", "B", "C"), n, replace = T)) 

# Simulate a phenotype composed of phylogenetic, fixed and residual effects 
d$phenotype <- rbv(bird.families, 1, nodes="TIPS") + 
       d$X1*0.7 + 
       ifelse(d$X3 == "B", 0.5, 0) + 
       ifelse(d$X3 == "C", 0.8, 0) + 
       rnorm(n, 0, 1) 

# Inverse matrix of shared phyloegnetic history 
Ainv <- inverseA(bird.families)$Ainv 

# Set priors 
prior <- list(R = list(V = 1, nu = 0.002), 
       G = list(G1 = list(V = 1, nu = 0.002))) 

# Run model 
model <- MCMCglmm(phenotype ~ X1 + X2 + X3 + X4, 
        random = ~taxon, 
        ginverse = list(taxon=Ainv), 
        data = d, 
        prior = prior, 
        verbose = FALSE) 

summary(model) 

# drop1(model, test = "Chisq") # doesn't work 
# update(model, ~.- X2)  # doesn't work 

model2 <- MCMCglmm(phenotype ~ X1 + X3 + X4, 
        random = ~taxon, 
        ginverse = list(taxon=Ainv), 
        data = d, 
        prior = prior, 
        verbose = FALSE) 

summary(model2) 

전체 모델 출력 :

Iterations = 3001:12991 
Thinning interval = 10 
Sample size = 1000 

DIC: 145.2228 

G-structure: ~taxon 

     post.mean l-95% CI u-95% CI eff.samp 
taxon  1.808 0.2167 3.347 17.88 

R-structure: ~units 

     post.mean l-95% CI u-95% CI eff.samp 
units 0.6712 0.0003382 1.585 20.98 

Location effects: phenotype ~ X1 + X2 + X3 + X4 

      post.mean l-95% CI u-95% CI eff.samp pMCMC  
(Intercept) -0.6005279 -1.5753606 0.2689091 1000.0 0.198  
X1   0.7506249 0.5220491 1.0033785 506.8 <0.001 *** 
X2   -0.0339347 -0.2452923 0.2085228 900.7 0.712  
X3B   0.6276852 0.1238884 1.1727409 1000.0 0.022 * 
X3C   1.0533919 0.4244092 1.5095804 1000.0 <0.001 *** 
X4B   -0.0002833 -0.5099753 0.5103896 478.3 0.978  
X4C   0.0723440 -0.4435520 0.6193011 1000.0 0.776  
--- 
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1 

감소 모델 출력 :

Iterations = 3001:12991 
Thinning interval = 10 
Sample size = 1000 

DIC: 296.4755 

G-structure: ~taxon 

     post.mean l-95% CI u-95% CI eff.samp 
taxon  1.542 0.07175  3.1 26.33 

R-structure: ~units 

     post.mean l-95% CI u-95% CI eff.samp 
units 0.8155 0.0004526 1.661 21.89 

Location effects: phenotype ~ X1 + X3 + X4 

      post.mean l-95% CI u-95% CI eff.samp pMCMC  
(Intercept) -0.606920 -1.380267 0.335225 1000.0 0.172  
X1   0.750929 0.540667 1.011041 1000.0 <0.001 *** 
X3B   0.633123 0.149282 1.243790 1000.0 0.032 * 
X3C   1.047085 0.518664 1.570641 1000.0 <0.001 *** 
X4B   0.007019 -0.528080 0.512161 509.2 0.998  
X4C   0.063093 -0.490103 0.590965 1000.0 0.818  
--- 
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1 

답변

2

모델 선정 과정에 상관없이, 매우 표준입니다 무엇을 당신이 적합하거나 사용하는 패키지의 종류.

다양한 후보 공변량 관련된 모든 회귀 모델

, 후진 선택이 재귀하지 :

1. fit a model with `p` covariates; 
2. drop the least insignificant covariate based on p-value (in your case, the pMCMC value); 
3. `p = p - 1` and go back to 1. 

The process goes on until there is no insignificant covariates to drop. 

마찬가지로 전방 절편으로부터 시작하여 연속적 공변량 추가 선택이있다. 나는 여기서 그것을 설명하지 않을 것이다.

모델 선택은 AIC/BIC와 같은 정보 기준을 기반으로 할 수도 있습니다. MCMC 추론의 경우 DIC가 사용됩니다. 가장 좋은 모델은 DIC가 가장 작은 모델입니다.

제가 제안하는 것은 이전 버전을 선택하는 것입니다. 다양한 모델을 피팅하는 동안 DIC 값도 기록하십시오. 결국 당신은 다음과 같은 테이블을 얻게됩니다 :

Model_1 DIC_1 
Model_2 DIC_2 
    .   . 
    .   . 

가장 좋은 상황은 DIC 선택이 동일한 모델을 제공한다는 것입니다. 이것은이 모델이 나머지 부분에 대한 강력한 증거를 가지고 있음을 의미합니다. 글쎄, 당신은 황금!

잔차, 원본 데이터 및 예상 표준 오류로부터 R- 제곱을 계산 (조정) 할 수도 있습니다. 가장 큰 R-square 모델이 가장 좋습니다.

+1

고마워요! 답은 단순화 프로세스에 접근하는 방법에 대한 개요를 제공하지만 update() 및 drop1()과 같은 함수가 없으면 몇 가지 예측 변수 이상으로 많은 노력이 필요할 것입니다. 따라서 간단한 접근 방법이 있는지 궁금해했습니다. 실제 코딩의 – user2390246

+1

물론 루프를 작성하는 것이 해결책이지만 사소한 것은 아닙니다. 특히 상호 작용이 시작되면 더욱 그렇습니다. 그러므로 더 쉬운 방법이 있는지 물어 봅니다! 그러나 예, 그럴 수 있습니다. 그리고 물론, 나는 그것을 누군가에게 나에게 쓰기를 요구하지는 않는다. – user2390246

+0

glmnet에 대해서는 잘 모르지만, MCMCglmm과 같이 계통 발생 정보를 포함 할 수 없다는 것을 알 수있다. – user2390246