R在水文时间序列分析的应用

R在水文时间序列分析的应用自回归滑动平均模型AutoregressiveModels-AR(p)ar{stats}FitAutoregressiveModelstoTimeSeriesDescriptionFitanautoregressivetimeseriesmodeltothedata,bydefaultselectingthecomplexitybyAIC.Usagear(x,aic=TRUE,order.max=NULL,method=c(yule-walker,burg,ols,mle,yw),na.action,series,...)ArgumentsxAunivariateormultivariatetimeseries.aicLogicalflag.IfTRUEthentheAkaikeInformationCriterionisusedtochoosetheorderoftheautoregressivemodel.IfFALSE,themodeloforderorder.maxisfitted.order.maxMaximumorder(ororder)ofmodeltofit.DefaultstothesmallerofN-1and10*log10(N)whereNisthenumberofobservationsexceptformethod=mlewhereitistheminimumofthisquantityand12.methodCharacterstringgivingthemethodusedtofitthemodel.Mustbeoneofthestringsinthedefaultargument(thefirstfewcharactersaresufficient).Defaultstoyule-walker.na.actionfunctiontobecalledtohandlemissingvalues.seriesnamesfortheseries.Defaultstodeparse(substitute(x)).在概率论中，一个时间序列是一串随机变量。在统计学中，这样一些变量都会受时间影响：比如每天在变的股票价格，每月一测的空气温度，每分钟病人的心率等等数据：北美五大湖之一的LakeHuron的1875-1972年每年的水位值这个时间序列大致的图像：plot(LakeHuron,ylab=,main=LevelofLakeHuron)AR（1）模型：x-LakeHuronop-par(mfrow=c(2,1))y-filter(x,.8,method=recursive)plot(y,main=AR(1),ylab=)acf(y,main=paste(p=,signif(dwtest(y~1)$p.value,3)))par(op)ACF和PCF图op-par(mfrow=c(3,1),mar=c(2,4,1,2)+.1)acf(x,xlab=)pacf(x,xlab=)spectrum(x,xlab=,main=)par(op)AR(p)模型使用Yule-walker法得出估计的参数值y-ar(x,aic=TRUE,method=yule-walker)regr=ar.ols(x,order=2,demean=FALSE,intercept=FALSE)regr结果：Call:ar.ols(x=x,order.max=2,demean=FALSE,intercept=FALSE)Coefficients:121.1319-0.1319Orderselected2sigma^2estimatedas0.5281预测1973值1.1319*x[98]-0.1319*x[97][1]579.9692参考书目：IntroductoryTimeSerieswithR，AnalysisofTimeSeriesDataUsingR，TimeSeriesAnalysisandItsApplications--withRexamples，TimeSeriesAnalysisandItsApplications--withRexamples参考网站：(MovingAveragemodels)Hereisasimplewayofbuildingatimeseriesfromawhitenoise:justperformaMovingAverage(MA)ofthisnoise.n-200x-rnorm(n)y-(x[2:n]+x[2:n-1])/2op-par(mfrow=c(3,1),mar=c(2,4,2,2)+.1)plot(ts(x),xlab=,ylab=whitenoise)plot(ts(y),xlab=,ylab=MA(1))acf(y,main=)par(op)n-200x-rnorm(n)y-(x[1:(n-3)]+x[2:(n-2)]+x[3:(n-1)]+x[4:n])/4op-par(mfrow=c(3,1),mar=c(2,4,2,2)+.1)plot(ts(x),xlab=,ylab=whitenoise)plot(ts(y),xlab=,ylab=MA(3))acf(y,main=)par(op)Youcanalsocomputethemovingaveragewithdifferentcoefficients.n-200x-rnorm(n)y-x[2:n]-x[1:(n-1)]op-par(mfrow=c(3,1),mar=c(2,4,2,2)+.1)plot(ts(x),xlab=,ylab=whitenoise)plot(ts(y),xlab=,ylab=momentum(1))acf(y,main=)par(op)n-200x-rnorm(n)y-x[3:n]-2*x[2:(n-1)]+x[1:(n-2)]op-par(mfrow=c(3,1),mar=c(2,4,2,2)+.1)plot(ts(x),xlab=,ylab=whitenoise)plot(ts(y),xlab=,ylab=Momentum(2))acf(y,main=)par(op)Insteadofcomputingthemovingaveragebyhand,youcanusethefilterfunction.n-200x-rnorm(n)y-filter(x,c(1,-2,1))op-par(mfrow=c(3,1),mar=c(2,4,2,2)+.1)plot(ts(x),xlab=,ylab=Whitenoise)plot(ts(y),xlab=,ylab=Momentum(2))acf(y,na.action=na.pass,main=)par(op)TODO:theside=1argument.AR(Auto-Regressivemodels)Anothermeansofbuildingatimeseriesistocomputeeachtermbyaddingnoisetotheprecedingterm:thisiscalledarandomwalk.Forinstance,n-200x-rep(0,n)for(iin2:n){x[i]-x[i-1]+rnorm(1)}Thiscanbewritten,moresimply,withthecumsumfunction.n-200x-rnorm(n)y-cumsum(x)op-par(mfrow=c(3,1),mar=c(2,4,2,2)+.1)plot(ts(x),xlab=,ylab=)plot(ts(y),xlab=,ylab=AR(1))acf(y,main=)par(op)Moregenerally,onecanconsiderX(n+1)=aX(n)+noise.Thisiscalledanauto-regressivemodel,orAR(1),becauseonecanestimatethecoefficientsbyperformingaregressionofxagainstlag(x,1).n-200a-.7x-rep(0,n)for(iin2:n){x[i]-a*x[i-1]+rnorm(1)}y-x[-1]x-x[-n]r-lm(y~x-1)plot(y~x)abline(r,col='red')abline(0,.7,lty=2)Moregenerally,anAR(q)processisaprocessinwhicheachtermisalinearcombinationoftheqprecedingtermsandawhitenoise(withfixedcoefficients).n-200x-rep(0,n)for(iin4:n){x[i]-.3*x[i-1]-.7*x[i-2]+.5*x[i-3]+rnorm(1)}op-par(mfrow=c(3,1),mar=c(2,4,2,2)+.1)plot(ts(x),xlab=,ylab=AR(3))acf(x,main=,xlab=)pacf(x,main=,xlab=)par(op)Youcanalsosimulatethosemodelswiththearima.simfunction.n-200x-arima.sim(list(ar=c(.3,-.7,.5)),n)op-par(mfrow=c(3,1),mar=c(2,4,2,2)+.1)plot(ts(x),xlab=,ylab=AR(3))acf(x,xlab=,main=)pacf(x,xlab=,main=)par(op)PACFThepartialAutoCorrelationFunction(PACF)providesanestimationofthecoefficientsofanAR(infinity)model:wehavealreadyseenitonthepreviousexamples.ItcanbeeasilycomputedfromtheautocorrelationfunctionwiththeYule-Walkerequations.Yule-WalkerEquationsTocomputetheauto-correlationfunctionofanAR(p)processwhosecoefficientsareknown,(1-a1B-a2B^2-...-apB^p)Y=Zwejusthavetocomputethefirstautocorrelationsr1,r2,...,rp,andthenusetheYule-Walkerequations:r(j)=a1r(j-1)+a2r(j-2)+...+apr(j-p).YoucanalsousethemintheotherdirectiontocomputethecoefficientsofanARprocessfromitsautocorrelations.DescriptionFitanARMAmodeltoaunivariatetimeseriesbyconditionalleastsquares.Forexactmaximumlikelihoodestimationseearima0.Usagearma(x,order=c(1,1),lag=NULL,coef=NULL,include.intercept=TRUE,series=NULL,qr.tol=1e-07,...)Argumentsxanumericvectorortimeseries.orderatwodimensionalinte