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）
（本资料素材和资料部分来自网络，仅供参考。
请预览后才下载，期待您的好评与关注！
）
升级会员 