前辈您好!为什么我运行你的程序出问题了呢?麻烦您帮我看下,谢谢您了!
> library(tsDyn)
载入需要的程辑包:mgcv
This is mgcv 1.7-3. For overview type 'help("mgcv-package")'.
载入需要的程辑包:Matrix
载入需要的程辑包:lattice
载入程辑包:'Matrix'
The following object(s) are masked from 'package:base':
det
载入需要的程辑包:snow
载入程辑包:'snow'
The following object(s) are masked from 'package:base':
enquote
载入需要的程辑包:mnormt
载入需要的程辑包:foreach
载入需要的程辑包:iterators
载入需要的程辑包:codetools
foreach: simple, scalable parallel programming from REvolution Computing
Use REvolution R for scalability, fault tolerance and more.
http://www.revolution-computing.com
载入需要的程辑包:MASS
> svpdx =read.table("quartc.txt", header = TRUE);
> x =svpdx$qcpi
> mod =star(x, m=5,d=1,thDelay=1,noRegimes=2,control=list(maxit=3000))
Testing linearity... p-Value = 0.002852162
The series is nonlinear. Incremental building procedure:
Building a 2 regime STAR.
Performing grid search for starting values...
Starting values fixed: gamma = 24 , th = 14.00072 ; SSE = 88.80465
Optimization algorithm converged
Optimized values fixed for regime 2 : gamma = 23.99999 , th = 13.99887
Finished building a MRSTAR with 2 regimes
>
> phi1=mod$model.specific$phi1
> phi1
[,1] [,2] [,3] [,4] [,5] [,6]
[1,] 0.1883161 1.4341176 -0.3197857 -0.13549745 -0.1913054 0.16428765
[2,] 14.4251639 -0.7906438 0.7040126 0.01523312 -0.5758762 -0.08938476
> phi2=mod$model.specific$phi2
> phi2
[,1] [,2]
[1,] 23.99999 13.99887
> re=mod$residuals
> n=length(x) #88
> n1=length(re) #83
> number=0*seq(1:n1)
> B=300
> gamma=0*seq(1:B)
>
> #Transition function
> G <- function(y, g, th) plogis(y, th, 1/g)
>
> for (b in 1:B) {
+ #bootstrap by resampling the residuals
+ #sampling with replacement from the original sample
+ rand=sample(seq(1:n1),n1, replace = TRUE)
+ rr=re[rand]
+ #generate new y
+ y=c(x[1:5], 6 :n)
+ for (t in 6:n){
+
+ y[t]=t(c(1,y[(t-1):(t-5)]))%*%phi1[1,] + (t(c(1,y[(t-1):(t-5)]))%*%phi1[2,])*G(y[t-2], phi2[1],phi2[2]) + rr[t-5]}
+
+ #{y[t]=0.18831613+1.43411762*y[t-1]-0.31978565*y[t-2]-0.13549745*y[t-3]-0.19130543*y[t-4]+0.16428765*y[t-5]+
+ #(14.42516390-0.79064384*y[t-1]+0.70401255*y[t-2]+0.01523312*y[t-3]-0.57587624*y[t-4]-
+ #0.08938476*y[t-5])/(1+exp(-23.99999222*(y[t-2]-13.99886588)))+re[t-5]} #generate new y
+ z=y
+ rs=star(z, m=5,d=1,thDelay=1,noRegimes=2,trace=FALSE,control=list(maxit=3000))
+ if (is.null(rs$coefficients[13]) == TRUE) gamma=NA else
+ gamma=rs$coefficients[13]
+ }
lstar: missing value during computations
lstar: missing value during computations
lstar: missing value during computations
lstar: missing value during computations
lstar: missing value during computations
警告信息:
1: In plogis(q, location, scale, lower.tail, log.p) : 产生了NaNs
2: In plogis(q, location, scale, lower.tail, log.p) : 产生了NaNs
3: In plogis(q, location, scale, lower.tail, log.p) : 产生了NaNs
4: In plogis(q, location, scale, lower.tail, log.p) : 产生了NaNs
5: In plogis(q, location, scale, lower.tail, log.p) : 产生了NaNs
> gamma=gamma[!is.na(gamma)] # remove NAs
> length(gamma) # 66
[1] 0
> length(gamma) # 66
[1] 0
> qq=quantile(gamma,probs=seq(0,1,0.05))
> qq
0% 5% 10% 15% 20% 25% 30% 35% 40% 45% 50% 55% 60% 65% 70% 75% 80% 85% 90% 95% 100%
NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA
63# epoh