1

我使用以下代码为每组观察值拟合不同的指数曲线，并且效果很好。 

p = c(10,20,15,25,20,30,25,35,30,40,25,35,20,30,15,25,10,20)
v = c(92,110,104,117,123,139,146,162,165,176,160,176,143,163,118,137,92,110)
group = factor(rep((1:9), each=2))

mm = model.matrix(~ 0 + group)

fit = nls(v ~ drop(mm %*% c(b1, b2, b3, b4, b5, b6, b7, b8, b9))*(1-exp(-k*p)),
      start = list(k=0.5, b1=1000, b2=2000, b3=3000, b4=4000, b5=5000, b6=6000, b7=7000, b8=8000, b9=9000))

summary(fit)
Formula: v ~ drop(mm %*% c(b1, b2, b3, b4, b5, b6, b7, b8, b9)) * (1 - 
exp(-k * p))

Parameters:
Estimate Std. Error t value Pr(>|t|)    
k    0.10928    0.01374   7.954 4.55e-05 ***
b1 129.13042    8.01108  16.119 2.20e-07 ***
b2 126.81086    6.43352  19.711 4.57e-08 ***
b3 141.74666    5.62817  25.185 6.61e-09 ***
b4 161.10250    5.06762  31.791 1.04e-09 ***
b5 174.94417    4.63884  37.713 2.68e-10 ***
b6 175.73100    5.20655  33.752 6.48e-10 ***
b7 165.58007    6.01256  27.539 3.26e-09 ***
b8 146.48962    6.92337  21.159 2.62e-08 ***
b9 129.13042    8.01108  16.119 2.20e-07 ***
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Residual standard error: 5.767 on 8 degrees of freedom

Number of iterations to convergence: 8 
Achieved convergence tolerance: 5.907e-06

对于不同的实验集（在 for 循环中），我需要多次运行此代码。棘手的部分是，对于每个实验，预测变量（即 b1、b2 等）的数量都会发生变化。

有没有简单的出路？

4

2 回答 2

0

您可以使用as.formula从字符串构建公式。

start <- as.list(c(k=0.5, setNames(1:9*1000, paste0("b", 1:9))))  # start values
groups <- list(g1=1:4, g2=2:4, g3=c(1,6:9))                       # some random groupings

lapply(groups, function(group){                                   # apply model to groups
    mm <- mm[, group]
    groupStr <- sprintf("b%d", group)
    nls(
        as.formula(paste("v~drop(mm%*%c(", paste(groupStr, collapse=","), "))*(1-exp(-k*p))")),
        start=start[c("k", groupStr)]
    )
})
于 2015-07-12T18:26:35.593 回答
0

您有一个部分线性的最小二乘模型。使用适当的算法，您不仅不必为您的组指定参数，而且还避免了对起始值的需要。

fit1 <- nls(v ~ mm * (1 - exp(-k * p)), algorithm = "plinear",
            start = list(k=0.5))
summary(fit1)

#Formula: v ~ mm * (1 - exp(-k * p))
#
#Parameters:
#             Estimate Std. Error t value Pr(>|t|)    
#k             0.10928    0.01374   7.954 4.55e-05 ***
#.lin.group1 129.13036    8.01107  16.119 2.20e-07 ***
#.lin.group2 126.81082    6.43351  19.711 4.57e-08 ***
#.lin.group3 141.74664    5.62816  25.185 6.61e-09 ***
#.lin.group4 161.10248    5.06761  31.791 1.04e-09 ***
#.lin.group5 174.94416    4.63883  37.713 2.68e-10 ***
#.lin.group6 175.73098    5.20655  33.752 6.48e-10 ***
#.lin.group7 165.58004    6.01256  27.539 3.26e-09 ***
#.lin.group8 146.48959    6.92337  21.159 2.62e-08 ***
#.lin.group9 129.13036    8.01107  16.119 2.20e-07 ***
#---
#Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
#
#Residual standard error: 5.767 on 8 degrees of freedom
#
#Number of iterations to convergence: 9 
#Achieved convergence tolerance: 6.958e-06
于 2015-07-12T19:43:55.117 回答