The code is cumbersome but looks right to me. However, I get "overflow" error when the program reaches about the 160th observation. I don't know if it's something wrong in my code or the simulated series is simply explosive. I doubt it's the latter since I can get the series simulated in R without problem. Thanks for your help!
Code: Select all
'y[i] <- rho1 + rho2*y[i-1] + rho3*y[i-2] + rho4*y[i-3] + rho5*y[i-4]+rho6*y[i-5] + rho7*y[i-6]+ rho8*y[i-7] + (rho9 + rho10*y[i-1] + rho11*y[i-2]+ rho12*y[i-3] + rho13*y[i-4] + rho14*y[i-5]+ rho15*y[i-6] + rho16*y[i-7])/(1+exp(-rho17*(y[i-1]-rho18)/0.32588125424785))+v[i-8]
‘blag_{%new} is the number of lags in the model, which is 7 in this example
‘The total observation is 182
series simulated_star_se=0 ‘generate an empty series as a starting point for simulation
for !i=1 to 182-blag_{%new}
!y=blag_{%new}+!i 'simulate current y, starting from 8th observation in this example to the end
for !p=1 to blag_{%new}
!y_lag=blag_{%new}+!i-!p 'lag in AR to simulate current y
!coef2=!p+2+blag_{%new} ‘This is c(10) when !p=1. Since c(1) is the constant in the 1st regime, c(9) is the constant in the 2nd regime, this is the coefficient for y(-1) in the 2nd regime
simulated_star_se({!y})=simulated_star_se({!y})+simulated_star_se({!y_lag})*coef_se({!coef2}) 'creating (rho9 + rho10*y[i-1] + rho11*y[i-2] + rho12*y[i-3] + rho13*y[i-4] + rho14*y[i-5]+ rho15*y[i-6] + rho16*y[i-7])
‘coef_se is the vector storing the coefficients from the STAR estimation
next
!rho17=2*{blag_{%new}}+3
!rho18=2*{blag_{%new}}+4
!lag_thresh=!y-1 ‘y(-1) is the threshold variable in this model, which is used in the transition equation
simulated_star_se({!y})=simulated_star_se({!y})/(1+exp(-coef_ex({!rho17})*(simulated_star_se({!lag_thresh})-coef_se(!rho18))/0.32588125424785)) 'times the above term by the transition function
for !m=1 to blag_{%new}
!y_lag=blag_{%new}+!i-!m
!coef1=!m+1
simulated_star_se({!y})=simulated_star_se({!y})+simulated_star_se({!y_lag})*coef_se({!coef1}) 'Add rho1 + rho2*y[i-1] + rho3*y[i-2] + rho4*y[i-3] + rho5*y[i-4]+rho6*y[i-5] + rho7*y[i-6] + rho8*y[i-7] term.
next
!constant1=1
!constant2=2+blag_{%new}
simulated_star_se({!y})=simulated_star_se({!y})+star_se_resid({!y})+coef_se({!constant1})+coef_se({!constant2})/(1+exp(-coef_se({!rho17})*(simulated_star_se({!lag_thresh})-coef_se({!rho18}))/0.32588125424785)) 'Add two constants and the residual term
next