R에서 nls를 사용하여 연구를 다시 작성하십시오.

Question

R에서 nls 함수를 사용하는 방법을 배우고 있으며 일부 문제가 있습니다. 나는 지금 연구 논문에서 발견 된 곡선을 재현하려고 시도하고있다. 이 모델은 1987 년 이전의 주식 시장 움직임에 곡선을 매 깁니다. 이 같은 NLS 불렀다

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

: 다음

I는 FUNC 함수를 정의

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

DATA2 두 열로 구성되는 데이터 프레임이고, 하나는이고 날짜가 다른 하나는 값입니다. 1466 개의 행이 있습니다. 내가 NLS를 실행할 때 다음과 같은 메시지를 얻을

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 

, 내가 무엇을 수집 할 수 있습니다에서

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 

,이 데이터 프레임이 설정되는 방식에 문제가 있지만 해결책을 찾을 수 없습니다.

내가 어떻게이 아버지를 따라갈 수 있을지 아십니까?

도움 주셔서 감사합니다.

Answer 1

기본 문제는 독립 변수를 지정하지 않았다는 것입니다. a, b, tc, and t에 start(...)을 지정하면 nls(...)에 모델의 모든 매개 변수가 표시됩니다.

a, b, and tc이 매개 변수이고 t이 독립 변수 인 LPPL 모델의 단순화 된 버전을 사용하는 것 같습니다. 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 = 1로 시작하는 반면에 b가 <으로 제한된다는 것을 언급합니다.] 그래서 어떻게해야합니까?

minpack 패키지의 minpack.lm(...) 기능은 비선형 최소 제곱 추정을 위해 매우 강력한 Levenberg-Marquardt 알고리즘을 사용합니다. 사실, 당신이 인용 한 논문은 L-M을 구체적으로 언급합니다. 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")) 

Answer 2

일반적으로는, 가정은 인 (귀하의 경우, Y) 변수 소위 "의존"또는 "응답"이 있다는 것입니다 하나 또는 그 이상의 "독립"또는 "예측 변수"(Date)의 함수이며 일반적으로 예측 함수 자체의 상세한 사양은 일반적으로 매우 적은 수의 정적 매개 변수 (a 및 b 및 아마도 t 및/또는 tc도 정확히 무엇인지에 따라 다시 달성하려고 노력한다). nls() 함수의 작업은 가장 정확한 가능한 예측을 유도하는 정적 매개 변수에 대한 최적 값을 찾는 것입니다.

예측 함수 func의 입력에 필요한 독립 변수가 누락 된 것 같습니다. 그래서 아마도 두 가지 중 하나를해야 할 것 같아요. func을 수정하여 Date을 입력으로 사용하거나 데이터 프레임의 Date 열 레이블을 변경하여 func 입력 중 하나와 일치하게하십시오 (대부분 Date 열의 이름을 바꿀 것으로 의심됩니다). 그것은 tc에 해당합니다). 두 경우 모두 고정 된 오프셋 날짜 (예 : (tc - t))에서 데이터 프레임의 날짜 값을 빼는 계산을 수행하려면 R이 실제로 있는지 확인해야합니다 날짜를 문자열이 아닌 Date 객체로 인식하여 의미있는 방식으로 다른 하나를 빼는 방법을 알고 있어야합니다. 이 목적으로 as.Date() 기능이 도움이 될 수 있습니다.

추가 대안으로, 오히려 입력으로 R 날짜 객체를 받아 너무 func를 다시 작성하는 것보다, 당신은 간단 단지를 참조하여 일의 경과 정수로 데이터 프레임에 Date 열을 재 할당 찾을 수 있습니다 약간의 상쇄; 예를 들면 다음과 같이하십시오 :

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

또는 그와 유사한 것.