5

考虑 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 或其他软件包之一中会更容易,那将是有用的知识,但是关于如何使用它们处理上述玩具问题的一些指示会有所帮助。

4

2 回答 2

9

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 的不同值自动执行此操作的方法。这使用文本修改构造公式,然后使用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))
}
于 2012-02-04T10:16:56.180 回答
2

我在澄清自己的问题时意识到,这让我看到有一种比我以前见过的更笨拙的方法。

即使有一些明显的流线型可以进入,这对我来说仍然有点不雅,但至少可以忍受重复使用,所以我认为它是一个足够的答案。我仍然对比下面这个更简洁的方式感兴趣。

Hong Ooi 在 ns 生成的矩阵上使用 data.frame 来自动命名列的技巧有点可爱,我在下面使用了它。一般来说,我可能会使用 paste 来构建它们,因为我有几个变量可以使用。

假设问题中给出的数据设置 -

lin.expr <- function(p,xn) {
  pn<-paste(p, 1:length(xn), sep = "")
  paste(paste(pn,xn,sep=" * "),collapse=" + ")
  }


m <- ns(x, df=3)
mydf <- data.frame(y, m)  # X-variables will be named X1, X2, ... 
xn <- names(mydf)[2:dim(mydf)[2]]

nspb <- lin.expr("b",xn)

c.form <- paste("y ~ theta * plogis( a + ",nspb,")",sep="")
stl <- list(theta=2, a=-5,b1=10, b2=10, b3=10)
nls( c.form, data=mydf, start= stl)

我的实际公式将包含几个术语,例如 nspb。实质性改进表示赞赏;我宁愿不选择自己的答案,但我想如果一两天内没有任何进一步的问题,我会选择它。

编辑:Hong Ooi 的补充(在我输入时发布并使用了类似的想法,但添加了一些不错的附加功能)几乎可以做到;这是一个可以接受的答案,所以我已经检查过了。

于 2012-02-05T03:17:22.163 回答