R中非线性最小二乘内的样条

Question

考虑R中的非线性最小二乘模型,例如以下形式:

 y ~ theta / ( 1 + exp( -( alpha + beta * x) ) )

(我的真正的问题有几个变量,外部函数不是逻辑但更多涉及;这个更简单,但我想如果我能做到这一点,我的情况应该几乎立即跟随)

我想用(例如)自然三次样条替换术语"alpha + beta*x".

这里有一些代码用于在逻辑内部创建一些非线性函数的示例数据:

set.seed(438572L)
x <- seq(1,10,by=.25)
y <- 8.6/(1+exp( -(-3+x/4.4+sqrt(x*1.1)*(1.-sin(1.+x/2.9))) )) + rnorm(x, s=0.2 )

如果我在lm中不需要逻辑,我可以轻松地用样条项替换线性项; 所以线性模型是这样的:

 lm( y ~ x ) 

然后成为

 library("splines")
 lm( y ~ ns( x, df = 5 ) )

生成拟合值很简单,并借助于(例如)rms包得到预测值似乎很简单.

实际上,将原始数据与基于lm的样条拟合拟合并不是太糟糕,但我有理由在逻辑函数中需要它(或者更确切地说,在我的问题中等效).

nls的问题是我需要为所有参数提供名称(我很高兴他们称之为(b1,...,b5)为一个样条拟合(并说c1,...,c6为另一个变量) - 我需要能够制作其中的几个).

是否有一种合理的方法来生成nls的相应公式,以便我可以用样条函数替换非线性函数内的线性项？

我能想到的唯一方法就是可以做到这一点有点尴尬和笨重,如果不编写一大堆代码就不能很好地概括.

(编辑以供澄清)对于这个小问题,我当然可以手工完成 - 写出由ns生成的矩阵中每个变量的内积的表达式,乘以参数的向量.但是,我必须为每个其他变量中的每个样条再次逐个编写整个项目,并且每次我更改任何样条曲线中的df时再次,并且如果我想使用cs而不是ns,则再次.然后,当我想尝试做一些预测(/插值)时,我们会得到一系列新的问题需要处理.我需要一遍又一遍地继续这样做,并且可能需要大量的结和几个变量,以便在分析后进行分析 - 我想知道是否有一种比写出每个单独术语更简洁,更简单的方法,无需编写大量代码.我可以看到一个相当牛逼的方式,这将涉及到相当多的代码,但是作为R,我怀疑有更简洁的方式(或更可能是3或4个更简洁的方式)只是躲避我.因此问题.

我以为我曾经看到有人在过去以相当不错的方式做过这样的事情,但对于我的生活,我现在找不到它; 我已经尝试了很多次来找到它.

[更具体地说,我通常希望能够尝试适合每个变量中的几个不同样条曲线 - 尝试几种可能性 - 以便看看我是否能找到一个简单的模型,但仍然适合这个目的是足够的(噪音真的非常低;合适的偏差可以达到很好的平滑效果,但只能达到一定程度).它更像是"找到一个漂亮的,可解释的,但足够的拟合函数",而不是任何接近推理和数据挖掘的东西都不是这个问题的真正问题.

或者,如果这比gnm或ASSIST或其他包装更容易,那将是有用的知识,但是关于如何继续上述玩具问题的一些指示将有所帮助.

Answer 1

ns实际上生成了一个预测变量矩阵.您可以做的是将该矩阵拆分为单个变量,并将它们提供给nls.

m <- ns(x, df=5)
df <- data.frame(y, m)  # X-variables will be named X1, ... X5
# starting values should be set as appropriate for your data
nls(y ~ theta * plogis(alpha + b1*X1 + b2*X2 + b3*X3 + b4*X4 + b5*X5), data=df,
        start=list(theta=1, alpha=0, b1=1, b2=1, b3=1, b4=1, b5=1))

ETA:这是针对不同df值自动执行此操作.这使用文本munging构造公式,然后用于do.call调用nls.警告:未经测试.

my.nls <- function(x, y, df)
{
    m <- ns(x, df=df)
    xn <- colnames(m)
    b <- paste("b", seq_along(xn), sep="")
    fm <- formula(paste("y ~ theta * plogis(1 + alpha + ", paste(b, xn, sep="*",
          collapse=" + "), ")", sep=""))
    start <- c(1, 1, rep(1, length=length(b)))
    names(start) <- c("theta", "alpha", b)
    do.call(nls, list(fm, data=data.frame(y, m), start=start))
}