ggplot:使用CI95%带拟合曲线(geom_smooth方法=“nls”)

4

这个this网站发布了有关酶动力学的数据:

Enz <- c("WT","WT","WT","WT","WT",
         "WT","WT","WT","WT","WT",
         "WT","WT","WT",
         "H297F","H297F","H297F",
         "H297F","H297F","H297F",
         "H297F","H297F")
S <- c(2.00, 1.00, 0.60, 0.50, 0.40, 
         0.30, 0.20, 0.10, 0.09, 0.08, 
         0.06, 0.04, 0.02, 
         0.05, 0.10, 0.20, 
         0.30, 0.40, 0.50, 
         1.00, 2.00)
v <- c(59.01, 58.29, 54.17, 51.82, 49.76, 
         45.15, 36.88, 26.10, 23.50, 22.26, 
         16.45, 13.67, 6.14, 
         11.8, 19.9, 30.3, 
         36.6, 40.2, 42.1, 
         47.8, 50.0)

绘图的代码:

ggplot(data=enzdata,         
            aes(x=S,            
            y=v,            
            colour = Enz)) +  
            geom_point() +            
            xlab("Substrate (mM)") +  
            ylab("Velocity (uM/min/mg.enzyme)") +    
            ggtitle("Glucose Dehydrogenase \n wild type and mutant") +  
            geom_smooth(method = "nls", 
                        formula = y ~ Vmax * x / (Km + x), 
                        start = list(Vmax = 50, Km = 0.2),
                        se = F, size = 0.5, 
                        data = subset(enzdata, Enz=="WT")) +
            geom_smooth(method = "nls", 
                        formula = y ~ Vmax * x / (Km + x), 
                        start = list(Vmax = 50, Km = 0.2),
                        se = F, size = 0.5, 
                        data = subset(enzdata, Enz=="H297F"))

enter image description here

您是否可以为两条独立曲线添加CI95%带?

使用@adiana的解决方案,它生成了下一个图:

enter image description here

1个回答

3

很不幸,由于nls预测有些棘手,无法使用ggplot2自动完成,您需要手动拟合模型。

首先拟合模型:

library(nls2)
nsmodel1<-nls(formula = v ~ Vmax * S / (Km + S),data=subset(enzdata, Enz=="WT"),start = list(Vmax = 50, Km = 0.2))
nsmodel2<-nls(formula = v ~ Vmax * S / (Km + S),data=subset(enzdata, Enz=="H297F"),start = list(Vmax = 50, Km = 0.2))

然后预测这两个区间。在这里找到as.lm.nls的代码。 http://www.leg.ufpr.br/~walmes/cursoR/ciaeear/as.lm.R
fit1<-predict(as.lm.nls(nsmodel1), interval = "confidence") 
fit2<-predict(as.lm.nls(nsmodel2), interval = "confidence") 
enzdata$lowerfit[enzdata$Enz=="WT"]<-fit1[,2]
enzdata$upperfit[enzdata$Enz=="WT"]<-fit1[,3]
enzdata$lowerfit[enzdata$Enz=="H297F"]<-fit2[,2]
enzdata$upperfit[enzdata$Enz=="H297F"]<-fit2[,3]

最后使用geom_ribbon绘制区间图,假设p是您先前的拟合结果。
p+geom_ribbon(aes(x=S,ymin=lowerfit,ymax=upperfit),data=subset(enzdata, Enz=="WT"),alpha=0.5)+
  geom_ribbon(aes(x=S,ymin=lowerfit,ymax=upperfit),data=subset(enzdata, Enz=="H297F"),alpha=0.5) 

as.lm.nls函数中,我不得不做出以下更改:fo <- as.formula(fo, env = proto:::as.proto.list(L)) - c0bra

网页内容由stack overflow 提供, 点击上面的
可以查看英文原文,
原文链接