2015-05-29 68 views
0

基本上我想编写一个程序,将随机化我的数据的顺序n次,然后完成一个生存分析并绘制输出过n编写一个程序来介绍排列

所以,让我们从采取以下通用数据matching()打包并创建治疗和未治疗人员的数据集。 Link to package

set.seed(123) 

library(Matching) 
data(lalonde) 

lalonde$age_cat <- with(lalonde, ifelse(age < 24, 1, 2)) 
attach(lalonde) 

lalonde$ID <- 1:length(lalonde$age) 


#The covariates we want to match on 
X = cbind(age_cat, educ, black, hisp, married, nodegr, u74, u75, re75, re74) 
#The covariates we want to obtain balance on 
BalanceMat <- cbind(age_cat, educ, black, hisp, married, nodegr, u74, u75, re75, re74, 
        I(re74*re75)) 
genout <- GenMatch(Tr=treat, X=X, BalanceMatrix=BalanceMat, estimand="ATE", M=1, 
        pop.size=16, max.generations=10, wait.generations=1) 
detach(lalonde) 

# now lets pair the the non-treated collisions to the treated 
# BUT lets pair WITHOUT REPLACEMENT 

mout <- Match(Y=NULL, Tr=lalonde$treat, X=X, 
       Weight.matrix=genout, M=2, 
       replace=FALSE, ties=TRUE) 

summary(mout) 
# we see that for 130 treated observations, we have 260 non-treated 
# this is because we set M=2 
# and yes length(lalonde$age[lalonde$treat==0]) == 260 but just follow me please 
# but this was done for a specific reason 

# now lets create a table for our 130+260 collisions 
treated <- lalonde[mout$index.treated,] 
# now we only want one occurence of the treated variables 
library(dplyr) 
treat_clean <- treated %>% 
    group_by(ID) %>% 
    slice(1) 

non.treated <- lalonde[mout$index.control,] 

# finally we can combine to form one clear data.set 
matched.data <- rbind(treat_clean, non.treated) 

现在,我们可以做一个条件Logistic回归确定或与re78(1987年赚来的钱)和治疗相关。为此我们需要生存包。 Link to package

library(survival) 

比方说,如果乘客收入在1978年

matched.data$success <- with(matched.data, ifelse(re78 > 8125, 1, 0)) 

output <- clogit(success ~ treat, matched.data, method = 'efron') 

summary(output) 

比8125更使我们看到,或用于治疗(治疗= 1)为1.495

发生了成功,我们可以保存为:

iteration.1 <- exp(output$coefficients[1]) 

现在我们从匹配包中读取(link)replace = FALSE请注意,如果FALSE, 匹配的顺序一般很重要。比赛将在 相同的顺序中找到的数据进行排序

所以我想要做的创建将用于n

  • 随机化拉隆德$ ID为了
  • 运行功能匹配处理
  • 运行clogit算法
  • 每次保存输出exp(output$coefficients[1])
  • 情节OR(exp(output$coefficients[1]))for each n

Essenece我想介绍排列到分析中。 如何才能做到这一点,当可以说N = 5

回答

1

我的replicate大风扇这样的事情:

X <- cbind(...)       # what you had before 
BalanceMat <- cbind(...)    # ditto 
lalonde$ID <- seq.int(nrow(lalonde)) 

results <- replicate(1000, { 
    ## not certain if it's just $ID order that matters 
    lalonde$ID <- sample(nrow(lalonde)) 
    ## lalonde <- lalonde[ sample(nrow(lalonde)), ] 

    ## ... 
    ## rest of your computation 
    ## ... 

    #### optionally return everything 
    ## output 
    #### return just the minimum 
    exp(output$coefficients[1]) 
}) 

#### if you returned output earlier, you'll need this, otherwise not 
## coef <- exp(sapply(results, function(z) z$coefficients[1])) 

## plot as needed 

我不知道你的意思只是ID事情的顺序,或者如果整个数据库的顺序;相应地调整replicate循环的第一对耦合线。

+0

\ replicate \很好,我喜欢你使用的方式\ sample(nrow(lalonde))\ – lukeg

1

您可以使用sample引进排列

data(lalonde) 
lalonde$age_cat <- with(lalonde, ifelse(age < 24, 1, 2)) 
lalonde$ID <- 1:length(lalonde$age) 
n <- 5 
res <- rep(NA, n) 
for (i in 1:n) { 
    lalonde <- lalonde[sample(1:nrow(lalonde)), ] # randomise order 
    ## rest of code 
    res[i] <- exp(output$coefficients[1]) 
} 

plot(1:n, res, main="Odds Ratios") 
+0

也是一个很好的答案,工作正常。感谢您的巨大贡献 – lukeg

+0

为什么在随机化订单时,您会使用'lalonde [样品(1:nrow(lalonde))']''而不是'lalonde [样品((lalonde))'] – lukeg

+0

感谢您的补充建议,如果你想从L1的'for'循环中省略'lalonde <-'(即写'lalonde [sample(1:nrow(lalonde)),]''而不是'lalonde < - lalonde [sample(1:nrow (lalonde)),]'随机化不会被应用,是正确的吗? – lukeg