jamehandsome

  •  
  • 2016年10月22日
  • 注册于 2014年1月24日
  • 直接install.packages("stringr"),再install.packages("ggplot2")

  • 回复 第2楼 的 doctorjxd:请问是以怎样的形式?箱式图?还是点连接?不好意思,我是新手哈。还望不吝赐教。谢谢

  • DID模型是在最小二乘法(Ordinary least squares,OLS)的基础上构建而成,其基本模型为:y=α+βt+γd+δtd+ε。求教在R中应该如何操作。谢谢,求高手啊!

  • High-Low Charts可以在SPSS中实现,请问能在R中实现吗?求代码……谢谢!

  • > library(R2WinBUGS)

    Loading required package: coda

    Loading required package: lattice

    Loading required package: boot

    Attaching package: ‘boot’

    The following object is masked from ‘package:lattice’:

    melanoma

    Attaching package: ‘R2WinBUGS’

    The following objects are masked from ‘package:R2OpenBUGS’:

    as.bugs.array, attach.all, attach.bugs, bugs, bugs.data, bugs.log,

    detach.all, detach.bugs, monitor, read.bugs, write.model

    > library(foreign)

    > REmodeldata<-read.dta("G:/normal document/statistics software and video/R/network meta/network/REmodel_data.dta")

    >

    > REmodel<-function(){

    + # Binomial likelihood, logit link

    + # Random effects model for multi-arm trials

    +

    + for(i in 1:ns){ # LOOP THROUGH STUDIES

    + w[i,1] <- 0 # adjustment for multi-arm trials is zero for control arm

    + delta[i,1] < - 0 # treatment effect is zero for control arm

    + mu ~ dnorm(0,.0001) # vague priors for all trial baselines

    + for (k in 1:na) { # LOOP THROUGH ARMS

    + r[i,k] ~ dbin(p[i,k],n[i,k]) # binomial likelihood

    + logit(p[i,k]) < - mu + delta[i,k] # model for linear predictor

    + rhat[i,k] <- p[i,k] * n[i,k] # expected value of the numerators

    +

    + #Deviance contribution

    + dev[i,k] < - 2 * (r[i,k] * (log(r[i,k]) -log(rhat[i,k]))

    + + (n[i,k]-r[i,k]) * (log(n[i,k]-r[i,k]) - log(n[i,k] -rhat[i,k]))) }

    +

    + # summed residual deviance contribution for this trial

    + resdev < - sum(dev[i,1:na])

    + for (k in 2:na) { # LOOP THROUGH ARMS

    + # trial-specific LOR distributions

    + delta[i,k] ~ dnorm(md[i,k],taud[i,k])

    + # mean of LOR di stributions (with multi-arm trial correction)

    + md[i,k] <- d[t[i,k]] - d[t[i,1]] + sw[i,k]

    + # precision of LOR distributions (with multi-arm trial correction)

    + taud[i,k] < - tau *2*(k -1)/k

    + # adjustment for multi-arm RCTs

    + w[i,k] <- (delta[i,k] - d[t[i,k]] + d[t[i,1]])

    + # cumulative adjustment for multi-arm trials

    + sw[i,k] < - sum(w[i,1:k -1])/(k-1)

    + }

    + }

    + totresdev < - sum(resdev[]) # Total Residual Deviance

    + d[1]<-0 # treatment effect is zero for reference t reatment

    + # vague priors for treatment effects

    + for (k in 2:nt){ d[k] ~ dnorm(0,.0001) }

    + sd ~ dunif(0,5) # vague prior for between -trial SD

    + tau <- pow(sd,-2) # between-trial precision = (1/between -trial variance)

    +

    + # pairwise ORs and LORs for all possi ble pair-wise comparisons, if nt>2

    + for (c in 1:(nt -1)) {

    + for (k in (c+1):nt) {

    + or[c,k] < - exp(d[k] - d[c])

    + lor[c,k] < - (d[k] -d[c])

    + }

    + }

    + # ranking on relative scale

    + for (k in 1:nt) {

    + rk[k] < - nt+1 -rank(d[],k) # assumes events are “good”

    + # rk[k] < - rank(d[],k) # assumes events are “bad”

    + best[k] < - equals(rk[k],1) #calculate probability that treat k is best

    + }

    + # *** PROGRAM ENDS

    + }# End of model

    >

    > filename <- file.path("D://","REmodel.bug")

    >

    > ## write model file:

    > write.model(REmodel, filename)

    > ## and let’s take a look:

    > file.show(filename)

    >

    > #load data from dataset

    >

    > t1<-REmodeldata$t1

    > r1<-REmodeldata$r1

    > n1<-REmodeldata$n1

    > t2<-REmodeldata$t2

    > r2<-REmodeldata$r2

    > n2<-REmodeldata$n2

    > t3<-REmodeldata$t3

    > r3<-REmodeldata$r3

    > n3<-REmodeldata$n3

    > na<-REmodeldata$na

    >

    > t<-c(t1,t2,t3)

    > r<-c(r1,r2,r3)

    > n<-c(n1,n2,n3)

    >

    > dim(t)<-c(25,3)

    > dim(r)<-c(25,3)

    > dim(n)<-c(25,3)

    > ns<-25

    > nt<-5

    >

    > data<-list("t","r","n","na","ns","nt")

    >

    > #Set Initial Values

    > inits<-function(){

    + #chain 1

    + list(d=c( NA, 0, 0, 0, 0), sd=1,

    + mu=c(0, 0, 0, 0, 0, 0, 0, 0, 0, 0,

    + 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,

    + 0, 0, 0, 0, 0))

    + #chain 2

    + list(d=c( NA, -1, -1, -1, -1) , sd=4,

    + mu=c(-3, -3, -3, -3, -3, -3, -3, -3, -3, -3,

    + -3, -3, -3, -3, -3, -3, -3, -3, -3, -3,

    + -3, -3, -3, -3, -3))

    +

    + #chain 3

    + list(d=c( NA, 2, 2, 2, 2), sd=2,

    + mu=c(-3, 5, -1, -3, 7, -3, -4, -3, -3, 0,

    + -3, -3, 0, 3, 5, -3, -3, -1, -3, -7,

    + -3, -3, 5, -1, 7))

    + }

    >

    > REmodel.sim<-bugs(data,inits,model.file="D:/REmodel.bug",

    + parameters=c("or","totresdev"),

    + n.chains=3,n.iter=40000,n.burnin=10000,

    + bugs.directory="D:/Program Files(x86)/WinBUGS14")

    Error in file(con, "rb") : cannot open the connection

    In addition: Warning message:

    In file(con, "rb") :

    cannot open file 'D:/Program Files(x86)/WinBUGS14/System/Rsrc/Registry.odc': No such file or directory

    Error in bugs.run(n.burnin, bugs.directory, WINE = WINE, useWINE = useWINE, :

    WinBUGS executable does not exist in D:/Program Files(x86)/WinBUGS14


  • > library(BRugs)

    > library(R2OpenBUGS)

    > library(foreign)

    > REmodeldata<-read.dta("G:/normal document/statistics software and video/R/network meta/network/REmodel_data.dta")

    > REmodel<-function(){

    + # Binomial likelihood, logit link

    + # Random effects model for multi-arm trials

    +

    + for(i in 1:ns){ # LOOP THROUGH STUDIES

    + w[i,1] <- 0 # adjustment for multi-arm trials is zero for control arm

    + delta[i,1] < - 0 # treatment effect is zero for control arm

    + mu ~ dnorm(0,.0001) # vague priors for all trial baselines

    + for (k in 1:na) { # LOOP THROUGH ARMS

    + r[i,k] ~ dbin(p[i,k],n[i,k]) # binomial likelihood

    + logit(p[i,k]) < - mu + delta[i,k] # model for linear predictor

    + rhat[i,k] <- p[i,k] * n[i,k] # expected value of the numerators

    +

    + #Deviance contribution

    + dev[i,k] < - 2 * (r[i,k] * (log(r[i,k]) -log(rhat[i,k]))

    + + (n[i,k]-r[i,k]) * (log(n[i,k]-r[i,k]) - log(n[i,k] -rhat[i,k]))) }

    +

    + # summed residual deviance contribution for this trial

    + resdev < - sum(dev[i,1:na])

    + for (k in 2:na) { # LOOP THROUGH ARMS

    + # trial-specific LOR distributions

    + delta[i,k] ~ dnorm(md[i,k],taud[i,k])

    + # mean of LOR di stributions (with multi-arm trial correction)

    + md[i,k] <- d[t[i,k]] - d[t[i,1]] + sw[i,k]

    + # precision of LOR distributions (with multi-arm trial correction)

    + taud[i,k] < - tau *2*(k -1)/k

    + # adjustment for multi-arm RCTs

    + w[i,k] <- (delta[i,k] - d[t[i,k]] + d[t[i,1]])

    + # cumulative adjustment for multi-arm trials

    + sw[i,k] < - sum(w[i,1:k -1])/(k-1)

    + }

    + }

    + totresdev < - sum(resdev[]) # Total Residual Deviance

    + d[1]<-0 # treatment effect is zero for reference t reatment

    + # vague priors for treatment effects

    + for (k in 2:nt){ d[k] ~ dnorm(0,.0001) }

    + sd ~ dunif(0,5) # vague prior for between -trial SD

    + tau <- pow(sd,-2) # between-trial precision = (1/between -trial variance)

    +

    + # pairwise ORs and LORs for all possi ble pair-wise comparisons, if nt>2

    + for (c in 1:(nt -1)) {

    + for (k in (c+1):nt) {

    + or[c,k] < - exp(d[k] - d[c])

    + lor[c,k] < - (d[k] -d[c])

    + }

    + }

    + # ranking on relative scale

    + for (k in 1:nt) {

    + rk[k] < - nt+1 -rank(d[],k) # assumes events are “good”

    + # rk[k] < - rank(d[],k) # assumes events are “bad”

    + best[k] < - equals(rk[k],1) #calculate probability that treat k is best

    + }

    + # *** PROGRAM ENDS

    + }# End of model

    >

    > filename <- file.path("D://","REmodel.bug")

    >

    > ## write model file:

    > write.model(REmodel, filename)

    > ## and let’s take a look:

    > file.show(filename)

    > #load data from dataset

    >

    > t1<-REmodeldata$t1

    > r1<-REmodeldata$r1

    > n1<-REmodeldata$n1

    > t2<-REmodeldata$t2

    > r2<-REmodeldata$r2

    > n2<-REmodeldata$n2

    > t3<-REmodeldata$t3

    > r3<-REmodeldata$r3

    > n3<-REmodeldata$n3

    > na<-REmodeldata$na

    >

    > t<-c(t1,t2,t3)

    > r<-c(r1,r2,r3)

    > n<-c(n1,n2,n3)

    >

    > dim(t)<-c(25,3)

    > dim(r)<-c(25,3)

    > dim(n)<-c(25,3)

    > ns<-25

    > nt<-5

    >

    > data<-list("t","r","n","na","ns","nt")

    >

    > #Set Initial Values

    > inits<-function(){

    + #chain 1

    + list(d=c( NA, 0, 0, 0, 0), sd=1,

    + mu=c(0, 0, 0, 0, 0, 0, 0, 0, 0, 0,

    + 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,

    + 0, 0, 0, 0, 0))

    + #chain 2

    + list(d=c( NA, -1, -1, -1, -1) , sd=4,

    + mu=c(-3, -3, -3, -3, -3, -3, -3, -3, -3, -3,

    + -3, -3, -3, -3, -3, -3, -3, -3, -3, -3,

    + -3, -3, -3, -3, -3))

    +

    + #chain 3

    + list(d=c( NA, 2, 2, 2, 2), sd=2,

    + mu=c(-3, 5, -1, -3, 7, -3, -4, -3, -3, 0,

    + -3, -3, 0, 3, 5, -3, -3, -1, -3, -7,

    + -3, -3, 5, -1, 7))

    + }

    >

    > REmodel.sim<-bugs(data,inits,model.file="D:/REmodel.bug",

    + parameters=c("or","totresdev"),

    + n.chains=3,n.iter=40000,n.burnin=10000,

    + bugs.directory="D:/Program Files(x86)/OpenBUGS/OpenBUGS323")

    Error in bugs(data, inits, model.file = "D:/REmodel.bug", parameters = c("or", :

    unused argument (bugs.directory = "D:/Program Files(x86)/OpenBUGS/OpenBUGS323")


  • 本人流行病与卫生统计小硕一枚,目前研二。正在学习R,希望大牛们不吝赐教!谢谢

  • 本人流行病与卫生统计学小硕一枚,目前研二。