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

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

《R在水文时间序列分析的应用.docx》由会员分享，可在线阅读，更多相关《R在水文时间序列分析的应用.docx（14页珍藏版）》请在冰豆网上搜索。

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

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

自回归滑动平均模型

AutoregressiveModels-AR（p）

ar{stats}FitAutoregressiveModelstoTimeSeries

Description

Fitanautoregressivetimeseriesmodeltothedata,bydefaultselectingthecomplexitybyAIC.

Usage

ar（x,aic=TRUE,order.max=NULL,

method=c（"yule-walker","burg","ols","mle","yw"）,

na.action,series,...）

Arguments

x

Aunivariateormultivariatetimeseries.

aic

Logicalflag.IfTRUEthentheAkaikeInformationCriterionisusedtochoosetheorderoftheautoregressivemodel.IfFALSE,themodeloforderorder.maxisfitted.

order.max

Maximumorder（ororder）ofmodeltofit.DefaultstothesmallerofN-1and10*log10（N）whereNisthenumberofobservationsexceptformethod="mle"whereitistheminimumofthisquantityand12.

method

Characterstringgivingthemethodusedtofitthemodel.Mustbeoneofthestringsinthedefaultargument（thefirstfewcharactersaresufficient）.Defaultsto"yule-walker".

na.action

functiontobecalledtohandlemissingvalues.

series

namesfortheseries.Defaultstodeparse（substitute（x））.

在概率论中，一个时间序列是一串随机变量。

在统计学中，这样一些变量都会受时间影响：

比如每天在变的股票价格，每月一测的空气温度，每分钟病人的心率等等

数据：

北美五大湖之一的LakeHuron的1875-1972年每年的水位值这个时间序列大致的图像：

plot（LakeHuron,

ylab="",

main="LevelofLakeHuron"）

AR

（1）模型：

x<-LakeHuron

op<-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:

12

1.1319-0.1319

Orderselected2sigma^2estimatedas0.5281

预测1973值

>1.1319*x[98]-0.1319*x[97]

[1]579.9692

参考书目：

IntroductoryTimeSerieswithR，AnalysisofTimeSeriesDataUsingR，TimeSeriesAnalysisandItsApplications--withRexamples，TimeSeriesAnalysisandItsApplications--withRexamples

参考网站：

http:

//zoonek2.free.fr/UNIX/48_R/15.html#2

MA（MovingAveragemodels）

Hereisasimplewayofbuildingatimeseriesfromawhitenoise:

justperformaMovingAverage（MA）ofthisnoise.

n<-200

x<-rnorm（n）

y<-（x[2:

n]+x[2:

n-1]）/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="MA

（1）"）

acf（y,main=""）

par（op）

n<-200

x<-rnorm（n）

y<-（x[1:

（n-3）]+x[2:

（n-2）]+x[3:

（n-1）]+x[4:

n]）/4

op<-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<-200

x<-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<-200

x<-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,youcanusethe"filter"function.

n<-200

x<-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:

the"side=1"argument.

AR（Auto-Regressivemodels）

Anothermeansofbuildingatimeseriesistocomputeeachtermbyaddingnoisetotheprecedingterm:

thisiscalledarandomwalk.

Forinstance,

n<-200

x<-rep（0,n）

for（iin2:

n）{

x[i]<-x[i-1]+rnorm

（1）

}

Thiscanbewritten,moresimply,withthe"cumsum"function.

n<-200

x<-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,onecanconsider

X（n+1）=aX（n）+noise.

Thisiscalledanauto-regressivemodel,orAR

（1）,becauseonecanestimatethecoefficientsbyperformingaregressionofxagainstlag（x,1）.

n<-200

a<-.7

x<-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<-200

x<-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）

Youcanalsosimulatethosemodelswiththe"arima.sim"function.

n<-200

x<-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）

PACF

ThepartialAutoCorrelationFunction（PACF）providesanestimationofthecoefficientsofanAR（infinity）model:

wehavealreadyseenitonthepreviousexamples.Itcanbeeasilycomputedfromtheautocorrelationfunctionwiththe"Yule-Walker"equations.

Yule-WalkerEquations

Tocomputetheauto-correlationfunctionofanAR（p）processwhosecoefficientsareknown,

（1-a1B-a2B^2-...-apB^p）Y=Z

wejusthavetocomputethefirstautocorrelationsr1,r2,...,rp,andthenusetheYule-Walkerequations:

r（j）=a1r（j-1）+a2r（j-2）+...+apr（j-p）.

YoucanalsousethemintheotherdirectiontocomputethecoefficientsofanARprocessfromitsautocorrelations.

Description

FitanARMAmodeltoaunivariatetimeseriesbyconditionalleastsquares.Forexactmaximumlikelihoodestimationsee arima0.

Usage

arma（x,order=c（1,1）,lag=NULL,coef=NULL,

include.intercept=TRUE,series=NULL,qr.tol=1e-07,...）

Arguments

x

anumericvectorortimeseries.

order

atwodimensionalintegervectorgivingtheordersofthemodeltofit. order[1] correspondstotheARpartand order[2] totheMApart.

lag

alistwithcomponents ar and ma.Eachcomponentisanintegervector,specifyingtheARandMAlagsthatareincludedinthemodel.Ifboth, order andlag,aregiven,onlythespecificationfrom lag isused.

coef

IfgiventhisnumericvectorisusedastheinitialestimateoftheARMAcoefficients.ThepreliminaryestimatorsuggestedinHannanandRissanen（1982）isusedforthedefaultinitialization.

include.intercept

Shouldthemodelcontainanintercept?

series

namefortheseries.Defaultsto deparse（substitute（x））.

qr.tol

the tol argumentfor qr whencomputingtheasymptoticstandarderrorsof coef.

Examples

data（tcm）

r<-diff（tcm10y）

summary（r.arma<-arma（r,order=c（1,0）））

summary（r.arma<-arma（r,order=c（2,0）））

summary（r.arma<-arma（r,order=c（0,1）））

summary（r.arma<-arma（r,order=c（0,2）））

summary（r.arma<-arma（r,order=c（1,1）））

plot（r.arma）

data（nino）

s<-nino3.4

summary（s.arma<-arma（s,order=c（20,0）））

summary（s.arma

<-arma（s,lag=list（ar=c（1,3,7,10,12,13,16,17,19）,ma=NULL）））

acf（residuals（s.arma）,na.action=na.remove）

pacf（residuals（s.arma）,na.action=na.remove）

summary（s.arma

<-arma（s,lag=list（ar=c（1,3,7,10,12,13,16,17,19）,ma=12）））

summary（s.arma

<-arma（s,lag=list（ar=c（1,3,7,10,12,13,16,17）,ma=12）））

plot（s.arma）

（本资料素材和资料部分来自网络，仅供参考。

请预览后才下载，期待您的好评与关注！

）