1

我觉得我已经接近找到我的问题的答案，但不知何故我无法做到。我已经使用 nls 函数来拟合 3 个参数，该函数使用一个相当复杂的函数来描述在一系列精子浓度（x 轴）中卵子的受精成功（y 轴）（Styan 模型[1]、[2]）。拟合参数效果很好，但我无法使用predict函数绘制平滑的外推曲线（见本文末尾）。我想这是因为我使用了一个未安装在 x 轴上的值。我的问题是如何根据nls在 x 轴上使用非拟合参数拟合函数的模型绘制平滑和外推曲线？

这是一个例子：

library(ggplot2)

data.nls <- structure(list(S0 = c(0.23298, 2.32984, 23.2984, 232.98399, 2329.83993, 
23298.39926), fert = c(0.111111111111111, 0.386792452830189, 
0.158415841584158, 0.898648648648649, 0.616, 0.186440677966102
), speed = c(0.035161615379406, 0.035161615379406, 0.035161615379406, 
0.035161615379406, 0.035161615379406, 0.035161615379406), E0 = c(6.86219803476946, 
6.86219803476946, 6.86219803476946, 6.86219803476946, 6.86219803476946, 
7.05624476582978), tau = c(1800, 1800, 1800, 1800, 1800, 1800
), B0 = c(0.000102758645352932, 0.000102758645352932, 0.000102758645352932, 
0.000102758645352932, 0.000102758645352932, 0.000102758645352932
)), .Names = c("S0", "fert", "speed", "E0", "tau", "B0"), row.names = c(NA, 
6L), class = "data.frame")

## Model S

modelS <- function(Fe, tb, Be) with (data.nls,{
x <- Fe*(S0/E0)*(1-exp(-B0*E0*tau))
b <- Fe*(S0/E0)*(1-exp(-B0*E0*tb))
x*exp(-x)+Be*(1-exp(-x)-(x*exp(-x)))*exp(-b)})

## Define starting values

start <- list(Fe = 0.2, tb = 0.1, Be = 0.1)

## Fit the model using nls

modelS.fitted <- nls(formula = fert ~ modelS(Fe, tb, Be), data = data.nls, start = start, 
control=nls.control(warnOnly=TRUE,minFactor=1e-5),trace = T, lower = c(0,0,0), 
upper = c(1, Inf, 1), algorithm = "port")

## Combine model parameters

model.data <- cbind(data.nls, data.frame(pred = predict(modelS.fitted)))

## Plot

ggplot(model.data) +
geom_point(aes(x = S0, y = fert), size = 2) +
geom_line(aes(x = S0, y = pred), lwd = 1.3) +
scale_x_log10()

在此处输入图像描述

我在这里尝试了 joran 的示例，但它没有效果，可能是因为我不适合S0：

r <- range(model.data$S0)
S0.ext <- seq(r[1],r[2],length.out = 200)
predict(modelS.fitted, newdata = list(S0 = S0.ext))
# [1] 0.002871585 0.028289057 0.244399948 0.806316161 0.705116868 0.147974213
4

1 回答 1

4

你的函数应该有参数(S0,E0,B0,tau,Fe,tb,Be)。nls将在传递给其data参数的 data.frame 中查找参数，并仅尝试适合那些在其中找不到的参数（前提是给出了起始值）。with在你的职能中不需要这个有趣的事情。（with无论如何都不应该在函数内部使用。它用于交互使用。） Inpredict newdata必须包含所有变量，即 S0、E0、B0 和 tau。

试试这个：

modelS <- function(S0,E0,B0,tau,Fe, tb, Be) {
  x <- Fe*(S0/E0)*(1-exp(-B0*E0*tau))
  b <- Fe*(S0/E0)*(1-exp(-B0*E0*tb))
  x*exp(-x)+Be*(1-exp(-x)-(x*exp(-x)))*exp(-b)}

## Define starting values

start <- list(Fe = 0.2, tb = 0.1, Be = 0.1)

## Fit the model using nls

modelS.fitted <- nls(formula = fert ~ modelS(S0,E0,B0,tau,Fe, tb, Be), data = data.nls, start = start, 
                     control=nls.control(warnOnly=TRUE,minFactor=1e-5),trace = T, lower = c(0,0,0), 
                     upper = c(1, Inf, 1), algorithm = "port")

## Combine model parameters

model.data <- data.frame(
  S0=seq(min(data.nls$S0),max(data.nls$S0),length.out=1e5),
  E0=seq(min(data.nls$E0),max(data.nls$E0),length.out=1e5),
  B0=seq(min(data.nls$B0),max(data.nls$B0),length.out=1e5),
  tau=seq(min(data.nls$tau),max(data.nls$tau),length.out=1e5))
model.data$pred <- predict(modelS.fitted,newdata=model.data)


## Plot

ggplot(data.nls) +
  geom_point(aes(x = S0, y = fert), size = 2) +
  geom_line(data=model.data,aes(x = S0, y = pred), lwd = 1.3) +
  scale_x_log10()

在此处输入图像描述

显然，这可能不是您想要的，因为该函数有多个变量并且不止一个变量在new.data. 通常，对于这样的情节，一个人只会改变一个，而其他人则保持不变。

所以这可能更合适：

  S0 <- seq(min(data.nls$S0),max(data.nls$S0),length.out=1e4)
  E0 <- seq(1,20,length.out=20)
  B0 <- unique(data.nls$B0)
  tau <- unique(data.nls$tau)

model.data <- expand.grid(S0,E0,B0,tau)
names(model.data) <- c("S0","E0","B0","tau")

model.data$pred <- predict(modelS.fitted,newdata=model.data)



## Plot

ggplot(model.data) +
  geom_line(data=,aes(x = S0, y = pred, color=interaction(E0,B0,tau)), lwd = 1.3) +
  geom_point(data=data.nls,aes(x = S0, y = fert), size = 2) +
  scale_x_log10()

在此处输入图像描述

于 2013-04-06T16:12:12.687 回答