2

我正在学习如何在 R 中使用 nls 函数并且遇到了一些问题。我现在只是试图重新创建一篇研究论文中发现的曲线。该模型拟合了 1987 年崩盘前股市走势的曲线。

我定义了一个函数func，如下：

func <- function(a,b,tc,t){
 a+b*log(tc-t)
}

我这样称呼 nls ：

nls1 <- nls(Y ~ func(a,b,tc,t), data2, start=list(a=0, b=1, tc=1466, t=1))

data2 是一个数据框，由两列组成，一列是日期，另一列是值。有 1466 行。

head(data2)
 Date      Y
1  1/4/82 882.52
2  1/5/82 865.30
3  1/6/82 861.02
4  1/7/82 861.78
5  1/8/82 866.53
6 1/11/82 850.46

我在运行 nls 时收到以下消息，

Error in qr(.swts * attr(rhs, "gradient")) : 
  dims [product 4] do not match the length of object [1466]

In addition: Warning message:

In .swts * attr(rhs, "gradient") :
  longer object length is not a multiple of shorter object length

据我所知，这是数据框设置方式的问题，但我找不到解决方案。

知道我怎样才能让这个父亲继续前进吗？

非常感谢您的帮助。

4

2 回答 2

8

基本问题是您没有指定自变量。通过指定start(...)for a, b, tc, and t，您是在告诉nls(...)这些都是模型的所有参数。

看起来您正在使用 LPPL 模型的简化版本，其中a, b, and tc是参数，并且t是自变量。它看起来像data2$Date包含时间变量。您需要确保data2$Date是 POSIXct 类。所以你可以写：

df$Date <- as.POSIXct(df$Date, format="%m/%d/%y")
nls1 <- nls(Y~a+b*log(tc-Date), data=data2, start=list(a=0, b=1, tc=1466))

编辑：回应OP的评论

这是一个很好的问题，因为它说明了使用nls(...). 您遇到的问题（现在已经正确指定了模型）nls(...)是没有收敛 - 一个令人痛苦的普遍现象。基本上，除非您的起始参数估计值相对接近最终的拟合值（或者除非模型非常“表现良好”），否则 nls 将失败。[还请注意，您引用的论文提到 b 限制为 b < 0，而您从 b = 1 开始。] 那该怎么办？

包中的minpack.lm(...)函数minpack使用异常稳健的 Levenberg-Marquardt 算法进行非线性最小二乘估计。事实上，你引用的论文特别提到了 LM。问题minpack.lm(...)在于它更难使用（您必须定义一个在给定步骤返回残差的函数，而不仅仅是定义要拟合的函数）。另外，minpack.lm(...)不计算拟合的统计数据。

所以解决方案是同时使用它们！用于minpack.lm(...)估计参数，然后将其用作 中的“起始值” nls(...)。下面的代码就是这样做的。使用拟合模型nls(...)将使生成拟合、预测值、残差的统计数据以及将模型应用于新数据集变得更加容易。

# this section just grabs the DJIA for 1982 - 1987; you already have this
library(tseries)
library(zoo)
ts <- get.hist.quote(instrument="DJIA", 
                     start="1982-01-01", end="1987-08-01", 
                     quote="Close", provider="yahoo", origin="1970-01-01",
                     compression="d", retclass="zoo")
df <- data.frame(ts)
df <- data.frame(Date=as.Date(rownames(df)),Y=df$Close)
df <- df[!is.na(df$Y),]
# end of setup...
library(minpack.lm) # for nls.lm(...)
library(ggplot2)    # for ggplot
df$days <- as.numeric(df$Date - df[1,]$Date)
# model based on a list of parameters
f <- function(pars, xx) {pars$a + pars$b*log(pars$tc - xx)} 
# residual function
resids <- function(p, observed, xx) {df$Y - f(p,xx)}
# fit using Levenberg-Marquardt algorithm
nls.out <- nls.lm(par=list(a=1,b=-1,tc=5000), fn = resids, 
                  observed = df$Y, xx = df$days)
# use output of L-M algorithm as starting estimates in nls(...)
par <- nls.out$par
nls.final <- nls(Y~a+b*log(tc-days),data=df, 
                 start=c(a=par$a, b=par$b, tc=par$tc))
summary(nls.final)      # display statistics of the fit 
# append fitted values to df
df$pred <- predict(nls.final)
# plot the results
ggplot(df)+
  geom_line(aes(x=Date,y=Y),color="black")+
  geom_line(aes(x=Date,y=pred),color="blue",linetype=2)+
  labs(title="LPPL Model Applied to DJIA (1982 - 1987)",
       x="", y="DJIA (daily close)")+
  theme(plot.title=element_text(face="bold"))

于 2013-12-30T05:33:24.300 回答
1

通常，当执行最小二乘回归时，假设存在一个所谓的“依赖”或“响应”变量（Y在您的情况下），它是一个或多个“独立”或“预测”的函数变量（Date），通常预测函数本身的详细规范通常由相当少量的静态参数定义（a和b，可能还有t和/或tc同样，取决于你想要实现的具体目标） . 该nls()函数的工作是为那些可能导致最准确预测的静态参数找到最佳值。

您的预测函数的输入func似乎缺少所需的自变量。所以，我认为你可能需要做两件事之一。要么进行修改func以使其接受Date作为输入，要么更改Date数据框中的列标签以使名称与输入之一匹配func（很可能我怀疑您希望重命名该Date列以使其对应于tc） . 在任何一种情况下，如果要执行计算，从固定偏移日期中减去数据框中的日期值（例如，(tc - t)就像现在写的那样），您将需要检查 R 是否实际上将您的日期识别为 Date 对象而不是字符串，以便它知道如何有意义地从另一个中减去一个。该as.Date()功能可能对您有所帮助。

作为另一种选择，与其尝试重写func以使其接受 R Date 对象作为输入，您可能会发现仅将Date数据框中的列重新分配给参考某个偏移量的经过的整数天数会更简单；例如，执行以下操作：

data2$tc <- as.numeric(as.Date(data2$Date) - as.Date("1982-1-4"))

或类似的。

于 2013-12-30T06:26:19.140 回答