不好意思,本人在统计和R语言上没什么经验,没办法简单的描述这个问题,如果各位觉得问题太浅,给个关键词我自己搜也行。谢谢。
我的目的
这是一个关于生物学活性计算的问题。我想用R语言复现文献中的SAS代码(原文如下)
proc nlin data=assay;
model obs=D+(A-D)/(1+(iu/((cs∗(sample=“S”)
+Ct∗(sample=“T”))))∗∗(B));
/∗语句中“∗”为乘法运算符,“∗∗”为指数运算符∗/
parms D=1 B=1 Cs=1 Ct=1 A=1;
run;
/∗assay数据集包含全部标准品和供试品的数据,
有3个变量,iu为浓度,obs为吸光度的观测值,
sample变量有两个取值,“S”和“T”分别表示该数据属
于标准品和供试品,
A、B、C、D分别为4个参数∗/
数据是这样的
dput(csf_1)
structure(list(conc = c(200, 100, 50, 25, 12.5, 6.25, 3.125,
1.5625, 200, 100, 50, 25, 12.5, 6.25, 3.125, 1.5625, 200, 100,
50, 25, 12.5, 6.25, 3.125, 1.5625, 200, 100, 50, 25, 12.5, 6.25,
3.125, 1.5625), id = c("s", "s", "s", "s", "s", "s", "s", "s",
"s", "s", "s", "s", "s", "s", "s", "s", "t", "t", "t", "t", "t",
"t", "t", "t", "t", "t", "t", "t", "t", "t", "t", "t"), abs = c(1.42,
1.408, 1.202, 0.84, 0.562, 0.423, 0.335, 0.312, 1.37, 1.338,
1.185, 0.843, 0.56, 0.391, 0.333, 0.302, 1.425, 1.395, 1.22,
0.863, 0.577, 0.413, 0.345, 0.317, 1.415, 1.364, 1.197, 0.862,
0.557, 0.404, 0.343, 0.313)), row.names = c(NA, -32L), class = "data.frame")
文献给出的参数值分别为:0.316(A)、1.833(B)、26.08(Cs)、25.11(Ct)、1.448(D)。我对s和t分别拟合也可以得出相近的数值。
我的代码和输出
经阅读《Nonlinear Regression with R》,我理解这个问题为固定a、b、d三个参数,c参数根据id变化,于是参考书中代码
Puromycin.m7 <- nls(rate ~ Vm[state] * conc/(K + conc), data = Puromycin, start = list(K = 0.1, Vm = c(200, 200)))
撰写以下代码
nlm_1 <- nls(abs ~ (a - d) / (1 + (conc / c[id]) ^ b) + d, data = csf_1, start = list(a = 0.3, b = 1.8, c = c(26, 26), d = 1.4))
但结果一直是
Error in numericDeriv(form[[3L]], names(ind), env) :
在计算模型的时候产生了缺省值或无限值
我不理解原因何在。请赐教。
系统环境
sessionInfo()
R version 4.0.0 (2020-04-24)
Platform: x86_64-w64-mingw32/x64 (64-bit)
Running under: Windows 10 x64 (build 17763)
Matrix products: default
locale:
[1] LC_COLLATE=Chinese (Simplified)_China.936
[2] LC_CTYPE=Chinese (Simplified)_China.936
[3] LC_MONETARY=Chinese (Simplified)_China.936
[4] LC_NUMERIC=C
[5] LC_TIME=Chinese (Simplified)_China.936
attached base packages:
[1] stats graphics grDevices utils
[5] datasets methods base
other attached packages:
[1] lattice_0.20-41 nlme_3.1-152 reshape2_1.4.4
[4] ggplot2_3.3.3
loaded via a namespace (and not attached):
[1] Rcpp_1.0.6 magrittr_2.0.1
[3] tidyselect_1.1.0 munsell_0.5.0
[5] colorspace_2.0-0 R6_2.5.0
[7] rlang_0.4.10 fansi_0.4.2
[9] stringr_1.4.0 plyr_1.8.6
[11] dplyr_1.0.5 tools_4.0.0
[13] grid_4.0.0 gtable_0.3.0
[15] utf8_1.2.1 DBI_1.1.1
[17] withr_2.4.1 ellipsis_0.3.1
[19] assertthat_0.2.1 tibble_3.1.0
[21] lifecycle_1.0.0 crayon_1.4.1
[23] purrr_0.3.4 vctrs_0.3.7
[25] glue_1.4.2 stringi_1.5.3
[27] compiler_4.0.0 pillar_1.5.1
[29] generics_0.1.0 scales_1.1.1
[31] jsonlite_1.7.2 pkgconfig_2.0.3