求助R软件中Logestic 回归分析的nomogram

是 否 +2 论坛币

k人 参与回答

经管之家送您一份

应届毕业生专属福利!

求职就业群

赵安豆老师微信：zhaoandou666

经管之家联合CDA

送您一个全额奖学金名额~ !

立即领取

感谢您参与论坛问题回答

经管之家送您两个论坛币！

+2 论坛币

       求助R软件中Logestic 回归分析的nomogram的制作和实施过程，因不是学统计的，所以R软件里的那个帮助内容没有看明白，请大侠给指个明路吧，十分感激～

扫码加我 拉你入群

请注明：姓名-公司-职位

以便审核进群资格，未注明则拒绝

分享0 收藏1 回帖

关键词：Nomogram 回归分析 gram OMO Log 回归分析 软件 统计 制作

话说如果没看懂，基本上是logistic回归需要回回炉～

  

这有一个例子，你看看能不能看懂——

  

n <- 1000    # define sample size

set.seed(17) # so can reproduce the results

age            <- rnorm(n, 50, 10)

blood.pressure <- rnorm(n, 120, 15)

cholesterol    <- rnorm(n, 200, 25)

sex            <- factor(sample(c('female','male'), n,TRUE))

# Specify population model for log odds that Y=1

L <- .4*(sex=='male') + .045*(age-50) +

  (log(cholesterol - 10)-5.2)*(-2*(sex=='female') + 2*(sex=='male'))

# Simulate binary y to have Prob(y=1) = 1/[1+exp(-L)]

y <- ifelse(runif(n) < plogis(L), 1, 0)

ddist <- datadist(age, blood.pressure, cholesterol, sex)

options(datadist='ddist')

f <- lrm(y ~ lsp(age,50)+sex*rcs(cholesterol,4)+blood.pressure)

nom <- nomogram(f, fun=function(x)1/(1+exp(-x)),  # or fun=plogis

    fun.at=c(.001,.01,.05,seq(.1,.9,by=.1),.95,.99,.999),

    funlabel="Risk of Death")

#Instead of fun.at, could have specified fun.lp.at=logit of

#sequence above - faster and slightly more accurate

plot(nom, xfrac=.45)

print(nom)

nom <- nomogram(f, age=seq(10,90,by=10))

plot(nom, xfrac=.45)

g <- lrm(y ~ sex + rcs(age,3)*rcs(cholesterol,3))

nom <- nomogram(g, interact=list(age=c(20,40,60)),

                conf.int=c(.7,.9,.95))

plot(nom, col.conf=c(1,.5,.2), naxes=7)

cens <- 15*runif(n)

h <- .02*exp(.04*(age-50)+.8*(sex=='Female'))

d.time <- -log(runif(n))/h

death <- ifelse(d.time <= cens,1,0)

d.time <- pmin(d.time, cens)

f <- psm(Surv(d.time,death) ~ sex*age, dist='lognormal')

med  <- Quantile(f)

surv <- Survival(f)  # This would also work if f was from cph

plot(nomogram(f, fun=function(x) med(lp=x), funlabel="Median Survival Time"))

nom <- nomogram(f, fun=list(function(x) surv(3, x),

                            function(x) surv(6, x)),

            funlabel=c("3-Month Survival Probability",

                       "6-month Survival Probability"))

plot(nom, xfrac=.7)

## Not run:

nom <- nomogram(fit.with.categorical.predictors, abbrev=TRUE, minlength=1)

nom$x1$points   # print points assigned to each level of x1 for its axis

#Add legend for abbreviations for category levels

abb <- attr(nom, 'info')$abbrev$treatment

legend(locator(1), abb$full, pch=paste(abb$abbrev,collapse=&rdquo;),

       ncol=2, bty='n')  # this only works for 1-letter abbreviations

#Or use the legend.nomabbrev function:

legend.nomabbrev(nom, 'treatment', locator(1), ncol=2, bty='n')

## End(Not run)

  

#Make a nomogram with axes predicting probabilities Y>=j for all j=1-3

#in an ordinal logistic model, where Y=0,1,2,3

Y <- ifelse(y==0, 0, sample(1:3, length(y), TRUE))

g <- lrm(Y ~ age+rcs(cholesterol,4)*sex)

fun2 <- function(x) plogis(x-g$coef[1]+g$coef[2])

fun3 <- function(x) plogis(x-g$coef[1]+g$coef[3])

f <- Newlabels(g, c(age='Age in Years'))  

#see Design.Misc, which also has Newlevels to change

#labels for levels of categorical variables

g <- nomogram(f, fun=list('Prob Y>=1'=plogis, 'Prob Y>=2'=fun2,

                     'Prob Y=3'=fun3),

         fun.at=c(.01,.05,seq(.1,.9,by=.1),.95,.99))

plot(g, lmgp=.2, cex.axis=.6)

options(datadist=NULL)

其他例子：http://biostat.mc.vanderbilt.edu/wiki/Main/DesignExamples