I'm trying to reproduce some research on the fitting LPPL on stock indices for bubble prediction and I'm having trouble with fitting the model to the data. I've been using the following papers for insight on this project: http://arxiv.org/pdf/1002.1010v2.pdf where they've already done some testing on the HSI, http://arxiv.org/pdf/0905.0220v1.pdf where I originally got my idea.
I've also tried to reproduce the findings from this stackoverflow post with no success (ran into multiple similar issues i.e. the max iterator has been reached, singular gradient matrix errors again): NLS And Log-Periodic Power Law (LPPL) in R
Not having much success with using daily prices to fit the model, I used weekly prices on the S&P following the advice in the conclusion of the HSI LPPL paper that the data should be "smoothed" in a way.
Here is the code I'm using. Advice on how to fix my issues would be much appreciated!
library(zoom)
library(minpack.lm)
library(tseries)
library(zoo)
#grab S&P500 historical
ts <- get.hist.quote(instrument="^GSPC",
start="2003-02-15", end="2007-10-31",
quote="Close", provider="yahoo", origin="1970-01-01",
compression="w", retclass="zoo")
df <- data.frame(ts)
df <- data.frame(Date=as.Date(rownames(df)),Y=df$Close)
df <- df[!is.na(df$Y),]
df$days <- as.numeric(df$Date - df[1,]$Date)
ts <- get.hist.quote(instrument="^GSPC",
start="1997-10-04", end="2011-10-12",
quote="Close", provider="yahoo", origin="1970-01-01",
compression="w", retclass="zoo")
df2 <- data.frame(ts)
df2 <- data.frame(Date=as.Date(rownames(df2)),Y=df2$Close)
df2 <- df2[!is.na(df2$Y),]
df2$days <- as.numeric(df2$Date - df2[1,]$Date)
f <- function(pars, xx)
with(pars,(a + ((tc - xx)^m) *b + c *(tc - xx)^m* cos(omega*log(tc - xx))+d *(tc - xx)^m* cos(omega*log(tc - xx))))
# residual function
resids <- function(p, observed, xx) {df$Y - f(p,xx)}
plot(df2$Date,df2$Y,type="l")
lines(df$Date,df$Y,type="l")
points(df$Date,df$Y,type="p")
pp = list(a=1662.239,b=-0.483332,tc=2050, m=0.97, omega=5, c=566, d=-566)
lines(df$Date,f(pars=pp,df$days),type="l")
nls.out <- nls.lm(par=pp, fn = resids, observed = df$Y, xx = df$days, control=nls.lm.control(maxiter=1000),lower = c(a = -Inf, b = -Inf, tc = 2008, m = 0.1, omega = 0.1, c = -Inf, d = -Inf), upper = c(a = Inf, b = -0.01, tc = 2050, m = 0.97, omega = 15, c = 3000, d = 3000))
par <- nls.out$par
par
lines(df$Date,f(par,df$days), col ="blue")
nls.out <- nls.lm(par=nls.out$par, fn = resids, observed = df$Y, xx = df$days, control=nls.lm.control(maxiter=1000),lower = c(a = -Inf, b = -Inf, tc = 2008, m = 0.1, omega = 3, c = -Inf, d = -Inf), upper = c(a = Inf, b = -0.01, tc = 2025, m = 0.999, omega = 10, c = Inf, d = Inf))
lines(df$Date,f(nls.out$par,df$days), col ="purple")
ppp = nls.out$par
lines(df$Date,f(ppp,df$days), col ="purple")
nls.final <- nls(Y~(a + ((tc - df$days)^m) * (b + c * cos(omega*log(tc - df$days))+d * cos(omega*log(tc - df$days)))), data=df, start=ppp, algorithm="port", control=nls.control(maxiter=1000, minFactor=1e-8), lower = c(a = -Inf, b = -Inf, tc = 2007, m = 0.01, omega = 6, c = -Inf, d = -Inf), upper = c(a = Inf, b = 0, tc = 2010, m = 0.999, omega = 10, c = Inf, d = Inf))
summary(nls.final) # display statistics of the fit
lines(df$Date,fitted(nls.final), col = "red")
# append fitted values to df
df$pred <- predict(nls.final, interval = "confidence")
summ = coef(summary(nls.final))
