R Time Series Analysis: Time Series Objects, ARIMA, and VAR

关于时间序列分析

时间序列分析的目的是帮助我们处理一类特殊的数据类型时间序列作为在经济金融领域建模最为常见的模型，在非常多的领域也有着应用；

时间序列分析中最为基础的手段自然是ARIMA模型以及从它衍生出的季节ARIMA模型；当然还有很多的用于不同的问题的模型后面一一介绍

时间序列数据类型与必要的基础

时间序列的数据可以保存在R的向量中，或者保存在R的数据框的一列或几列中，对应的时间单独保存或者保存在同一数据框中；

不过由于时间序列数据的特殊性 R提供了一些专用的格式来保存时间序列数据，他们本质上是一种特殊的数据框，但是提供了更方便交互的语法

ts类别是R基本的时间序列类别很多函数都基于ts形式实现在需要调用的时候我们应该用函数 as.ts 讲数据的格式进行转换 zoo形式是ts形式的一个扩展，允许的不等间隔的实现序列；整个zoo形式被xts形式吸收，更多高级的时间序列分析函数都基于xts形式构建 我们一般以xts格式存储数据，在需要的时候转换为ts格式使用

ts类型

创建ts类型序列

ts是基本R软件的stats包中支持的规则时间序列类型，具有start和frequency两个属性

日数据最好不时间标签为具体日历时间的ts类型，因为金融数据的日数据在周末和节假日一般无数据，而ts类型要求时间是每一天都相连的，不能周五直接跳到周一。

当然如果我们直接连接让交易日数据直接连接也未尝不可，这样可能在多元方面会出现一些问题因为有些数据可能在非交易日也可以收集

ts类型用于保存一元或者多元的等间隔时间序列，如月度、季度、年度数据，生成方法例子为

1	ts(x, start=c(2001, 1), frequency=12)

其中x是向量或者矩阵，取矩阵值时矩阵的每一列是一个序列。 frequency对月度数据是12，对季度数据是4，对年度数据可以缺省（值为1）

ts类型一般不用于日观测数据的情况，如果一定要这样使用那就当作年度数据设置频率为1

ts序列的常用函数

用ts.intersect函数和ts.union函数可以把两个或多个时间序列合并为多元时间序列，时间段取交集或者取并集,

为了使用序列数据进行计算、绘图等，可以用as.vector把一元时间序列的数据转换成普通向量

对于多元时间序列x，可以用x[,1]这样的格式取出其中分量时间序列，可以用as.vector(x[,1])将其中的分量转换为普通向量，可以借助于xts类型用coredata(as.xts(x))将多元时间序列的数据转换为普通矩阵。

用start()求时间序列的开始点， end()求时间序列的结束点， frequency()求采样频率

aggregate()函数可以把月度数据加总成年数据（实际上是一种数据降频）他的作用和普通数据框的aggregate()函数作用类型按照某种标签计算某种统计量时间序列将这种分类变量设置为了时间

time()函数对ts类型数据返回序列中的每个时间点的时间，结果是一个和原来时间序列形状相同的时间序列。 cycle()函数对月度数据返回序列每个时间点所在的月份, 结果是和原序列时间点相同的一个时间序列

window()函数取出时间序列的一段，如果指定frequency=TRUE还可以仅取出某个月（季度）

filter函数可以计算递推的或卷积的滤波类似zoo类型与xts类型的apply系列函数

zoo类型

R的zoo扩展包提供了比基本R中ts时间序列类型更灵活的时间序列类型，其时间标签（time stamp）可以使用R中任何的日期和时间类型，序列不需要是等时间间隔的，支持多元时间序列。

如果序列符合ts类型的要求，则与ts类型兼容并可以互相转换。 zoo也尽可能提供与ts类相同或相似功能的函数。 zoo扩展包提供的时间序列数据类型叫做zoo类型

作为对原本ts形式的扩展，我们可以大胆的混用其中的语法规则与函数，直到出错了再单独进行考虑

创建zoo类型序列

生成zoo类型的时间序列，只要提供两部分输入：一个向量或者矩阵x作为观测值，一个下标序列order.by作为排序变量以及时间标签。 不接受数据框

zoo的设计特点是允许使用任何可排序的数据类型作为时间标签

1	zoo(x, order.by)

继承 ts 类型的思想我们每一行对应一个固定的时间节点这意味着一个多元时间序列共用时间不过我们可以接受NA的存在所以这种灵活性已经足够了

举出例子

## 一元的时间序列
set.seed(1)
z.1 <- zoo(sample(3:6, size=12, replace=TRUE),
           make_date(2018, 1, 1) + ddays(0:11));
## 多元的时间序列
set.seed(2)
z.2 <- zoo(cbind(x=sample(5:10, size=12, replace=TRUE),
                 y=sample(8:13, size=12, replace=TRUE)),
           make_date(2018, 1, 1) + ddays(0:11));

能看出我们的order.by部分用一些函数生成了时间标签后面会介绍他们

ts、irts（在tseries包中定义）、its等时间序列类型可以用as.zoo()转换为zoo类型（极强的包容性），如果zoo类型的时间下标符合转换要求，也可以将zoo时间序列用as.xxx()类函数转换为其它时间序列类型。

文本文件中、字符串中、数据框中保存的时间序列数据可以用read.zoo()转换为zoo时间序列。

zoo类型泛型函数扩展

zoo类型算是比较核心的一个时间序列扩展类型他扩展了一些函数

print(x)显示x，对一元序列横向显示，对多元序列纵向显示

str(x)显示x的类型和结构信息这是对原本的函数str的扩展

head(x)和tail(x)可以取出序列开头的若干项与末尾的若干项也是扩展泛型函数

summary(x)对每个序列以及时间下标作简单概括统计也是扩展泛型函数

zoo类型子集提取

我们说过时间序列的封装都是扩展数据框形式因此提取语法可以参考数据框与矩阵形式进行

无论是一元还是多元序列提取方式没有变换

## 提取第一行
z[1]
## 提取一些行
z[10:12]
#提取一些行和第二列
z.2[1:3, 2]
## 使用时间下标进行提取行
z.2[ymd(c("2018-01-01", "2018-01-12"))]

日期时间字符串的格式为CCYY-MM-DD HH:MM:SS，而且可以省略后面的一部分，其含义是取出前面部分能匹配的所有时间点

zooreg类型

为了与规则时间序列的ts类型相对应， zoo扩展包提供了zoo的派生类型zooreg，代表规则间隔的时间序列，它具有与ts相同的时间间隔信息，但允许内部某些时间点不存在

用zoo()函数加frequency选项或者 zooreg()函数生成zooreg类型的时间序列语言基本继承

1 2	zoo(x, order.by, freqency) zooreg(x, start, end, frequency, deltat, ts.eps, order.by)

zooreg类型比ts优点是允许非NA的缺失同时保持类型不变 同时默认支持了以天为单位的观测

可以用is.regular(x)判断一个zoo类型数据是否规则，用is.regular(x, strict=TRUE)判断一个zoo类型数据是否满足和ts一样条件的规则序列

用as.ts(x)将一个规则的zooreg或者zoo类型数据转换为ts类型的时间序列。如果x内部有缺失的时间点，转换为ts类型时将补上这些时间点，观测值填以NA 用as.zoo(x)将ts类型的时间序列转换为zoo类型

zoo类型函数

zoo类型合并

两个保存了不同时间段的时间序列x, y可以用c(x,y)或者 rbind(x,y)合并，两个不同属性的时间序列，可以合并为一个多元时间序列，时间取并集，不存在的值取缺失值。两个序列不需要等长

zoo类型降维

类似于ts时间序列，对zoo时间序列，也可以用aggregate()对数据进行降频把多个小时间的观测值整合到一个较长时间里对时间序列而言这属于频率的改变

1	z.apy <- aggregate(z.ap, year, sum)

也是对函数aggregate的扩展

zoo类型填补

当序列存在缺失的时候我们有简易的填补方法，用na.locf()可以填补缺失值，填前面最近的一个非缺失值，用na.approx()可以对缺失值填充为线性插值，更多的方法查文档获得

zoo类型修改

两个序列之间可以进行加减乘除四则运算和逻辑比较运算，结果是对应时间点的相应运算。缺失值被删除

zoo时间序列类型的实际数据可以用coredata()读取或修改。读取时对一元序列结果是数值型向量，对多元序列结果是矩阵我们使用对矩阵的操作就可以实现原始数据的修改如

1
2
3

z.2b <- z.2
coredata(z.2b) <- 100 + coredata(z.2b)
z.2b

zoo类型的滑动平均

对时间序列经常计算像七日均线这样的滑动平均。对于一般的向量，基本R的函数filter()可以进行各种加权滑动平均和自回归迭代计算。 zoo包对zoo时间序列和ts时间序列提供了rollapply()函数，可以计算包括滑动平均在内的多种滚动计算

对一些常用的滚动计算zoo包含提供了单独的函数，如 rollmean(), rollmedian(), rollmax(), rollsum()，这几个是结果时间点对应滚动窗口中心的。而rollmeanr(), rollmedianr(), rollmaxr(), rollsumr() 则是结果时间点对准滚动窗口末尾时间点的

xts类型

继承与发展

xts包的设计目标是使得以xts为输入的函数可以兼容其它时间序列类型的输入，并可以很方便地增加属性。事实上xts已经成为了时间序列分析的第一优先级

xts包提供了xts时间序列类型，这本质上是zoo包的zoo类型，所以针对zoo类型的做法也适用于xts类型。他只是在zoo类型上进行了很小的修改

xts类型的创建

可以用xts()生成新的xts时间序列类型的数据对象，用法类似用zoo::zoo()

其它时间序列类型可以用as.xts()转换成xts类型

xts类型数据的子集

他无缝的衔接zoo类型子集提取的所有方法

可以用"from/to"的格式指定一个日期时间范围，而且也不需要开始点和结束点恰好有数据如

1	xts.1["2018-01-10/2018-01-14"]

first(x, n)和last(x, n)类似于head(x, n)和last(x, n)，但是对xts对象x，n除了可以取正整数值以外，还允许用字符串指定时间长度，允许的单位包括 secs, seconds, mins, minutes, hours, days, weeks, months, quarters, years。比如

1	first(xts.ap, "3 months")

字符串中取负值时表示扣除

xts类型函数

xts类型扩展函数

扩展了泛型函数 plot 查询plot.xts即可

可以用coredata(x)返回x的不包含时间的纯数据；用index(x)返回x的时间标签

periodicity(x)求xts对象x的时间区间

endpoints(x, on)给出按某种频率的分界点，频率包括"us"(微秒), "microseconds", "ms"(毫秒), "milliseconds", "secs", "seconds", "mins", "minutes", "hours", "days", "weeks", "months", "years"

xts金融时间序列降维

对于OHLC形式的金融时间序列，即开盘价(open)、最高价(high)、最低价(low)、收盘价(close)四个成分的时间序列，并且成分变量名也采用Open, High, Low, Close，可以用to.period(x, period)将其降频成period指定的采样频率 根据金融时间序列的习惯，选择收盘价是默认的选择

如果x是分钟数据， to.minutes3(x), to.minutes5(x), to.minutes10(x), to.minutes15(x),to.minutes30(x), to.hourly(x) 将x降频为3分钟到60分钟的数据

to.daily(x)将x降频为每天的数据，并在时间下标中删去时间部分。 to.weekly(x)将x降频为每周的数据，并在时间下标中删去时间部分，日期采用每周最后一天（周六）的日期。to.monthly(x)将x降频到月度数据并将时间下标改为yearmon类型， toquarterly(x)将x降频到季度数据并将时间下标改为yearqtr类型。 toyearly(x)将x降频到年度数据，日期采用每年有数据的最后一天的日期

xts类型滑动平均

用period.apply(x, INDEX, FUN)可以用INDEX对时间序列x分组，对每组用FUN函数计算一个函数值。 INDEX常取为endpoints(x, on="...")，给出某种分组周期。比如

1
2
3

period.apply(xts.ap,
             INDEX=endpoints(xts.ap, on="years"),
             FUN=mean)

常用操作如求和等有专门的函数，如period.sum(), period.min(), period.max(), period.prod()。专用函数运行效率更高。

常用的周期也有专门的函数，如apply.daily(), apply.weekly(), apply.monthly(), apply.quarterly(), apply.yearly()。这些函数实际是调用period.apply()

quantmod包

quantmod包的目的是为量化投资者提供方便的原型开发测试工具，而不是提供新的统计方法。

quantmod包提供了金融时间序列数据的一些方便功能，比如从公开数据源载入数据，作股票行情图、时间序列图等

quantmod包提供了一个getSymbols()函数，可以从多个公开数据源下载金融和经济数据并转换为R格式（主要是xts格式）

chartSeries()可以做K线图和曲线图封装了各种专门用于分析金融时间序列的功能是金融时间序列的优秀可视化方法

线性时间序列模型

对于线性时间序列的理论部分我们推荐线性时间序列分析同时 R语言的研究实例文件夹中存储了一份关于线性时间序列分析的例子解释参考它可以让让我们理解整个线性时间序列分析的算法实现细节

基本处理

无论我们使用什么类型的时间序列对象（当然我们推荐使用xts类型，它兼容性极好）这些基本处理函数都是通用的

扩展泛型函数 plot() 做最基本的时间序列曲线

由 stats 包提供函数 acf 给出样本自相关图， forecast包提供了一个类似功能的Acf()函数它没有保留0阶滞后可以让我们专注于后面的相关性 pacf函数则是计算样本偏自相关图，forecast包提供Pacf有一样的作用这些函数也会给出自相关系数的值给我们访问

lag()可以计算滞后序列，对ts类型输入， lag()的作用是序列数值不变，但是时间标签增加一个单位或者用k=指定的间隔

1 2	## 这才是传统意义上的滞后一个单位需要设置-1 x2 <- stats::lag(x1, k=-1); x2

diff(x)计算一阶差分， diff(x, lag, differences)计算滞后为lag的阶数位differences的差分也就是季节差分需要考虑的问题

ARIMA杂项

arima.sim可以模拟生成ARIMA模型的数据

补充ARFIMA fracdiff::fdGPH()计算差分阶的Geweke-Porter-Hudak估计值 fracdiff::fracdiff()函数进行ARFIMA模型估计

ARIMA模型辨识

ARIMA模型的基本辨识需要使用ACF和PACF曲线辅助在基本处理这里复用前面基本处理中的要点：acf() 用来查看样本自相关图，pacf() 用来查看样本偏自相关图；forecast::Acf() 和 forecast::Pacf() 提供类似功能，并且更适合在建模辨识时快速查看相关结构。

stats包中的ar()函数可以对时间序列样本进行AR建模默认采用AIC确定模型阶数提供了多种方法用于参数估计

forecast包提供了一个auto.arima()函数，可以自动进行模型选择，不过经常被人诟病效果不好

TSA包提供的eacf()函数辨识模型采用EACF方法线性时间序列分析：ARMA模型的识别

TSA包还提供了armasubsets()函数用来选择ARMA模型阶

ARIMA模型拟合

arima()函数估计一般性的ARIMA模型但是它需要预先指定阶数默认采用MLE估计模型阶数 arima函数允许指定某些系数固定为预先确定的值使用参数fixed确定，我们一般会直接把那些比较小的系数设置为0 从而实现稀疏估计 arima()函数可以用seasonal=指定季节模型，包括季节AR阶、季节差分阶、季节MA阶以及周期 arima()函数提供了一个 xreg= 用来引入回归自变量，他的做法是进行原序列作为因变量引入的回归自变量作为自变量残差使用ARMA模型分析可以看出 arima()模型是非常常用的时间序列分析函数对他的作用整合有

支持平稳可逆ARMA建模，这时可以用include.mean=TRUE设定有均值参数。支持ARIMA建模，如果差分阶大于零，则不允许include.mean=TRUE，即单位根过程不允许带有漂移。支持带有季节项的ARIMA建模。季节差分阶大于零时也不允许有漂移。支持以平稳ARMA序列为误差项的回归建模。

在arima()函数中用xreg=引入回归自变量（外生变量）时，不要指定差分或者季节差分

forecast包的Arima()函数能实现stats::arima()的功能，但在差分阶大于零时允许有漂移项

ARIMA模型检测

对于残差的正态性分析方法在其他地方都介绍过了，这里不重复叙述

Box.test()执行Ljung-Box白噪声检验它检验序列是否是白噪声序列对于残差分析这是非常有用的因为残差是从模型估计计算得到的，自由度有损失用fitdf=指定自由度减少个数选择为 $p+q$

可以用forecast::checkresiduals()函数进行模型诊断它同时进行残差分析的所有常见分析

将arima()的输出结果输入到tsdiag()函数，可以进行模型诊断是直接对模型的系统性诊断

单位根检验是一种平稳性检验，零假设是有单位根，即不平稳；对立假设是平稳。经常使用增强的Dickey-Fuller检验(ADF检验)；fUnitRoots包的adfTest()函数可以执行单位根ADF检验。 tseries包的adf.test()函数也可以执行单位根ADF检验; 单位根检验选项type选择基础模型可以取：

"nc"，表示没有漂移项或截距项；
"c"，表示带有一个漂移项或截距项；
"ct"，表示基础模型中带有 $a+bt$ 这样的线性项;

ARIMA模型预测

模型的预测还是使用经典的泛型函数 predict() 这是 stats提供的最基本方法它会提供预测配套的SE

用forecast包的forecast 函数进行预测是更加方便的

时间序列分解

stats包的decompose()函数输入一个时间序列，将其分解为趋势项、季节项和随机项。去趋势的方法是中心对称滑动平均。可以用type="additive"或type="multiplicative" 指定各项之间是相加还是相乘

stats包提供了函数stl()，该函数基于loess局部加权回归估计季节项，可以减少异常点的影响，属于稳健回归。用同月份（季度）的数值估计平滑的季节变动，减去季节项后再用平滑方法估计趋势

stats包的StructTS()函数用状态空间模型表示时间序列分解，用最大似然方法估计各个成分

stats包的HoltWinters()函数提供指数平滑方法，forecast包的ets()函数提供了自动选择合适的指数平滑方法并进行预报的功能。

ARCH系列模型

ARCH效应的检验

想要检验ARCH效应我们有两种方法

对残差的平方进行白噪声检验使用函数 Boxtest()
使用最小二乘方法对残差进行检验有函数FinTS::ArchTest() 注意，第一个函数输入残差序列平方，第二个直接接受残差序列

我们也有直观的检验方法也就是研究ACF曲线；我们要求残差项本身是白噪声序列（检验或者ACF都能体现这一点）而残差项的平方有自相关性（检验或者ACF都能体现这一点）这就是完全的检验方法

ARCH系列模型拟合

fGarch包的garchFit()可以函数建立ARCH模型他是波动率建模中最常见的函数 fGarch包的garchFit()函数支持多种条件分布，默认为正态分布，用cond.dist=指定分布：

"norm": 正态；
"snorm": 有偏正态;
"ged": 广义误差分布；
"sged": 有偏广义误差分布;
"std": t分布;
"sstd": 有偏t分布；
"snig"
"QMLE"：拟最大似然估计，仍假设正态但是采用稳健标准误差估计；选定条件分布以后再进行正态性检验就没有意义了我们的目标也不是正态分布

更复杂的ARCH类型建模需要使用扩展包 rugarch

ARCH模型的检验

对模型的检测就是研究残差序列的特征 residuals(model, standardize=TRUE) 函数可以计算标准化残差 summary(model) 可以快捷的给出模型概括其中包括了所有常规检测 plot.garch() 扩展了泛型函数plot 可以绘制所有常规的检测图形

ARCH拟合与预测

volatility() 函数可以拟合波动率 fitted() 函数拟合模型本身的值也就是均值项和波动率扩展泛型函数 predict() 函数用于多步预测也是扩展泛型函数

两步估计法

从原理上就可以看出我们不需要使用除了线性时间序列模型以外的内容这里需要我们自行设计函数来实现需要的效果原理介绍参考金融时间序列分析(一元)：两步估计法

多元时间序列分析与协整分析

多元时间序列分析主要使用的R扩展包有 MTS vars

模型的估计

MTS包的VAR()函数用于估计VAR模型，其第二个参数可以设定模型的阶数

MTS包的VARorder函数可以计算VAR定阶的 $M(i)$ 统计量和各种信息准则

vars包的VAR()函数也可以用于VAR模型的估计，其表现形式和MTS包不同，但是计算的结果一致，这个函数也允许自动选阶数

MTS包的refVAR()函数输入无约束的VAR建模结果，以及thres=1.645或thres=1.96这样的 $t$ 比值界限，生成设置部分系数为零的约束估计结果

模型的检验

MTS包的mq()函数用于进行多元混成检验，也就是一元中的Ljung-Box检验，允许手动设置自由度扣除。

自由度扣除根据系数来确定，使用残差进行多元混成检验的时候，有 $k^2p$ 个系数被估计，也就是我们设定的自由度扣除数量；如果存在一些系数被简化为0，他们不被包含在自由度扣除的范围中

MTS包还提供了一个MTSdiag()函数，输入模型结果和adj=自由度缩减个数，作残差的CCM估计表（ACF变种）、图和残差的多元混成检验

MTS包中GrangerTest()函数执行格兰杰因果性检验，在GrangerTest()中，用locInput=输入一个分量序号，检验零假设：其它所有分量都不是此分量的格兰杰原因。不能进行其它的格兰杰原因的检验。此函数的使用文档写作并不友好

模型的预测

MTS包的VARpred()函数可以从VAR的建模结果计算点预测值，不考虑参数估计误差的预测标准误差（Standard errors of predictions），考虑参数估计误差的预测标准误差（Root mean squared errors of predictions）

协整分析与与向量修正模型

扩展包tseries中po.test()可以执行基于EG两阶段法步骤的Phillips-Ouliaris协整检验，零假设是非协整，对立假设是存在协整关系

扩展包urca的ca.jo()函数可以进行计算Johansen的两种检验

状态空间模型

R软件中有许多扩展包支持进行状态空间模型建模。

statespacer：支持线性高斯状态时间序列建模。
KFAS：支持线性高斯状态时间序列建模，也支持指数分布族的非高斯情况。
dlm：使用(West and Harrison 1997)的模型，动态线性模型。支持线性高斯时间序列模型，可使用最大似然估计，支持时变模型。
dynr：支持带有机制切换的离散时间或者连续时间的模型建模。
dse: 线性高斯的ARMA、VAR和状态空间模型，所用的方法不太符合R的使用习惯。
bssm: 非线性、非高斯状态空间模型的贝叶斯推断。
MARSS: 多元自回归状态空间模型。参数可以时变，观测值可以包含缺失值。
MSwM: 马尔可夫机制切换的一元自回归模型，支持线性和广义线性模型。我们这里仅仅介绍其中的一部分

statespacer

statespacer包支持线性高斯状态空间模型建模，文档比较详细，使用的模型记号与初始化思想可以参考金融时间序列分析(一元)：线性高斯状态空间模型

对常用的结构时间序列模型、ARMA等提供了更简单的设定功能。对各矩阵时变的情况支持不足。

statespacer()是主要的建模函数。对于常见模型，可以使用选项直接指定模型。需要输入超参数的初始值， 这需要熟悉模型如何用状态空间表示。

结果都存成了深入嵌套的列表，设sr为保存了statespacer()结果的列表，访问其中的成分如下：

system_matrices: 系统矩阵
- H, Z, T, R, Q等。
predicted: 一步预测分布
- yfit: 的一步预测
- v: 的一步预测的误差
- Fmat: 的一步预测的误差方差阵
- a: 的一步预测
- P: 的一步预测的误差方差阵
- …………
filtered: 的滤波分布均值a和方差阵P。
smoothed: 平滑结果
- a: 的平滑分布的均值
- V: 的平滑分布的方差阵
- …………

MARSS

参考的模型结构表示为金融时间序列分析(一元)：MARSS包的模型,MARSS扩展包有详细的数百页的用户手册，在需要的时候我们可以参考。

MARSS包的基本函数是MARSS()，用法为

1	fit <- MARSS(y, model=list(...))

其中 y 是要建模的时间序列，如果是一维的也可以是普通R向量，更一般地可以输入为一个 $n\times T$ 矩阵，其中 $n$ 是观测值 $\boldsymbol{y}_t$ 的维数， $T$ 是时间序列观测时间点个数，每一列对应于一个时间点。观测值允许有缺失值。

model 用一个列表输入 B , U, 0 等成分，变量名就用金融时间序列分析(一元)：MARSS包的模型中的记号，但 $\boldsymbol\pi$ 用 x0 表示。如果还指定了 V0 , 则表示 $\boldsymbol{x}_{0}$ 先验分布指定为均值 x0 ,方差阵 V0

MARSS()返回一个marseMLE类型的对象，可以用各种信息提取函数提取或者进一步分析：

print(MLEobj)显示主要结果。 summary(MLEobj)显示结果更少一些。
coef(MLEobj)提取参数估计。用tidy::broom(MLEobj)提取成数据框格式。
residuals(MLEobj)提取观测或者状态的预测、滤波或平滑的残差，返回数据框形式。
tsSmooth(MLEobj)提取预测、滤波或平滑结果，用type参数选择，默认为平滑，返回数据框形式。加选项interval = "confidence"可以同时输出平滑的预测区间。 fitted(MLEobj)默认为一步预测，并且在预测、滤波、平滑时不会对噪声部分进行估计，所以要拟合观测因变量值应该使用此函数。支持对缺失数据的估计。在预测时可以用n.ahead指定步数作多步预测。
logLik(MLEobj)返回对数似然函数值。
AIC(MLEobj)返回AIC值， AICc(MLEobj)是对小样本情形进行修正的一个变种。 MLEobj <- MARSSaic(MLEobj)则添加更多的AIC类判别准则值。
MARSSkf(MLEobj)进行滤波、平滑，结果中： xtt1是状态的一步预报期望； xtt是状态的滤波期望； xtT是状态的平滑期望。 Vtt1, Vtt, VtT分别是状态的一步预报、滤波、平滑的方差阵估计，等等，见MARSS用户手册§3.3和§5.10。
MARSSparamCIs()计算参数置信区间，默认使用海色阵计算，用method = "parametric"指定使用参数bootstrap方法，用method = "innovation"指定使用新息重抽样方法。对于方差阵参数，不应使用海色阵方法。 bootstrap方法计算时间很长，应仅在算例充足时使用。
MARSSboot(MLEobj)用bootstrap方法进行置信区间、偏差估计等，可以用参数方法，或者新息重抽样方法，当观测值包含缺失值时仅支持新息重抽样方法。
还有一些进一步计算的函数，见MARSS文档中用户手册的§2.4。

隐马氏模型

可用的R扩展包：

depmixS4;
HiddenMarkov;
msm;
R2OpenBUGS(用于贝叶斯估计);
HMM(仅支持类别值时间序列)。