澳门新浦京娱乐场网站-www.146.net-新浦京娱乐场官网
做最好的网站

澳门新浦京娱乐场网站时间序列,时间序列分析

目录

目录

目录

目录


  • 时间序列分析工具箱——tibbletime
    • tibbletime 的用途
    • 加载包
    • 数据
    • 教程:tibbletime
      • 初始化一个 tbl_time 对象
      • 时间序列函数
  • 时间序列分析工具箱——timetk
    • timetk 的主要用途
    • 加载包
    • 数据
    • timetk 教程:
      • PART 1:时间序列机器学习
      • PART 2:转换
  • 时间序列分析工具箱——sweep
    • sweep 的用途
    • 加载包
    • 数据
    • Demo:forecast sweep 的简化预测工作流
      • STEP 1:创建 ts 对象
      • STEP 2A:ARIMA 模型
      • STEP 2B:简化模型
      • STEP 3:预测
      • STEP 4:用 sweep 简化预测
      • STEP 5:比较真实值和预测值
  • 时间序列分析工具箱——tidyquant
    • tidyquant 的用途
    • 加载包
    • tq_get:获得数据
      • 从 Yahoo! Finance 获得股票数据
      • 从 FRED 获得经济数据
    • 使用 tq_transmutetq_mutate 转换数据
      • tq_transmute
      • tq_mutate
    • 可用函数

第1章 准备工作
第2章 Python语法基础,IPython和Jupyter
第3章 Python的数据结构、函数和文件
第4章 NumPy基础:数组和矢量计算
第5章 pandas入门
第6章 数据加载、存储与文件格式
第7章 数据清洗和准备
第8章 数据规整:聚合、合并和重塑
第9章 绘图和可视化
第10章 数据聚合与分组运算
第11章 时间序列
第12章 pandas高级应用
第13章 Python建模库介绍
第14章 数据分析案例
附录A NumPy高级应用
附录B 更多关于IPython的内容(完)

翻译自《Demo Week: Tidy Time Series Analysis with tibbletime》

原文链接:www.business-science.io/code-tools/2017/10/26/demo_week_tibbletime.html

注意:由于软件包的版本变化,部分代码被修改,文字有删减

翻译自《Demo Week: Time Series Machine Learning with timetk》

原文链接:www.business-science.io/code-tools/2017/10/24/demo_week_timetk.html

翻译自《Demo Week: Tidy Forecasting with sweep》

原文链接:www.business-science.io/code-tools/2017/10/25/demo_week_sweep.html

时间序列分析工具箱——tidyquant

本文翻译自《Demo Week: class <- tidyquant》

原文链接:

澳门新浦京娱乐场网站 1


时间序列分析工具箱——tibbletime

澳门新浦京娱乐场网站 2

时间序列分析工具箱——timetk

澳门新浦京娱乐场网站 3

时间序列分析工具箱——sweep

澳门新浦京娱乐场网站 4

tidyquant 的用途

使用 tidyquant 的六大理由:

  1. 直接从 Yahoo! Finance、FRED Database、Quandl 等数据源获得网络数据
  2. 简化 xts、zoo、quantmod、TTR 和 PerformanceAnalytics 中金融及时间序列函数的调用
  3. 可视化: 漂亮的主题以及针对金融的 geom(例如 geom_ma
  4. 构建投资组合
  5. 财务分析以及投资组合归因方法
  6. 为金融与时间序列分析提供坚实的基础tidyquant 会自动加载 tidyverse 和各种金融、时间序列分析包,这使得它成为任何金融或时间序列分析的理想起点。

该教程将会介绍前两个主题。其他主题请查看 tidyquant 的文档。

时间序列(time series)数据是一种重要的结构化数据形式,应用于多个领域,包括金融学、经济学、生态学、神经科学、物理学等。在多个时间点观察或测量到的任何事物都可以形成一段时间序列。很多时间序列是固定频率的,也就是说,数据点是根据某种规律定期出现的(比如每15秒、每5分钟、每月出现一次)。时间序列也可以是不定期的,没有固定的时间单位或单位之间的偏移量。时间序列数据的意义取决于具体的应用场景,主要有以下几种:

tibbletime 的用途

  1. tidy 时间序列分析的未来:基于 tbl 的新类——tbl_time,为 tibble 对象添加时间轴,赋予处理时间的能力。
  2. 时间序列函数:为 tbl_time 对象专门设计的一系列函数,例如:
    • filter_time():根据日期简便快捷地过滤一个 tbl_time 对象。
    • as_period():转换时间周期(例如月度变为年度),让用户能将数据聚合到低粒度水平上。
    • time_collapse():当使用 time_collapse 时,tbl_time 对象中落入相同周期的索引将被修改成相同的日期。
    • rollify():修改一个函数,使其能够在特定时间区间上计算一个或一组值。可以用来计算滚动均值,或其他 tidyverse 框架下的滚动计算。
    • create_series():根据规则时间序列,用简化标记快速初始化一个带有 datetbl_time 对象。

timetk 的主要用途

三个主要用途:

  1. 时间序列机器学习:使用回归算法进行预测;
  2. 构造时间序列索引:基于时间模式提取、探索和扩展时间序列索引;
  3. 转换不同类型的时间序列数据(例如 tblxtszoots 之间):轻松实现不同类型的时间序列数据之间的相互转换。

我们今天将讨论时间序列机器学习和数据类型转换。第二个议题(提取和构造未来时间序列)将在时间序列机器学习中涉及,因为这对预测准确性非常关键。

sweep 的用途

正如 broom 包之于 stats 包,sweep 包用来简化使用 forecast 包的工作流。本教程将逐一介绍常用函数 sw_tidysw_glancesw_augmentsw_sweep 的用法。

sweeptimetk 带来的额外好处是,如果 ts 对象是由 tbl 对象转换来的,那么在预测过程中日期和时间信息会以 timetk 索引的形式保留下来。一句话概括:这意味着我们最终可以在预测时使用日期,而不是 ts 类型数据使用的规则间隔数字日期。

加载包

请先安装 tidyquant

# Install librariesinstall.packages("tidyquant")

加载 tidyquant

# Load librarieslibrary(tidyquant) # Loads tidyverse, financial pkgs, used to get and manipulate data
  • 时间戳(timestamp),特定的时刻。
  • 固定时期(period),如2007年1月或2010年全年。
  • 时间间隔(interval),由起始和结束时间戳表示。时期(period)可以被看做间隔(interval)的特例。
  • 实验或过程时间,每个时间点都是相对于特定起始时间的一个度量。例如,从放入烤箱时起,每秒钟饼干的直径。

加载包

tibbletime 目前还在活跃开发阶段,可以用常规方法安装,也可以借助 devtools 从 github 上安装最新开发版。

# Get tibbletime version with latest features
devtools::install_github("business-science/tibbletime")

安装完成后,加载下面的包:

  • tibbletime:创建带时间轴的 tibble 对象,可以使用 tbl_time 函数。
  • tidyquant:加载 tidyverse 框架,用 tq_get() 获取数据。
# Load libraries
library(tibbletime) # Version: 0.1.1, Future of tidy time series analysis
library(tidyquant)  # Loads tidyverse, tq_get()

加载包

需要加在两个包:

  • tidyquant:用于获取数据并在后台加载 tidyverse
  • timetk:R 中用于处理时间序列的工具包

如果还没有安装过,请用下面的命令安装:

# Install packagesinstall.packagesinstall.packages("tidyquant")

加载包。

# Load librarieslibrary     # Toolkit for working with time series in Rlibrary(tidyquant)  # Loads tidyverse, financial pkgs, used to get data

加载包

本教程要使用到四个包:

  • sweep:简化 forecast 包的使用
  • forecast:提供 ARIMA、ETS 和其他流行的预测算法
  • tidyquant:获取数据并在后台加载 tidyverse 系列工具
  • timetk:时间序列数据处理工具,用来将 tbl 转换成 ts
# Load librarieslibrary      # Broom-style tidiers for the forecast packagelibrary   # Forecasting models and predictions packagelibrary(tidyquant)  # Loads tidyverse, financial pkgs, used to get datalibrary     # Functions working with time series

tq_get:获得数据

使用 tq_get() 获得网络数据。tidyquant 提供了大量 API 用于连接包括 Yahoo! Finance、FRED Economic Database、Quandl 等等在内的数据源。

本章主要讲解前3种时间序列。许多技术都可用于处理实验型时间序列,其索引可能是一个整数或浮点数(表示从实验开始算起已经过去的时间)。最简单也最常见的时间序列都是用时间戳进行索引的。

数据

tq_get() 下载 FANG(脸书、亚马逊、网飞、谷歌)每天的股票价格。

# Stock Prices from Yahoo! Finance
FANG_symbols <- c("FB", "AMZN", "NFLX", "GOOG")

FANG_tbl_d <- FANG_symbols %>%
    tq_get(
        get = "stock.prices",
        from = "2014-01-01",
        to = "2016-12-31") 

FANG_tbl_d <- FANG_tbl_d %>%
    group_by(symbol)

FANG_tbl_d

## # A tibble: 3,024 x 8
## # Groups:   symbol [4]
##    symbol       date  open  high   low close   volume adjusted
##     <chr>     <date> <dbl> <dbl> <dbl> <dbl>    <dbl>    <dbl>
##  1     FB 2014-01-02 54.83 55.22 54.19 54.71 43195500    54.71
##  2     FB 2014-01-03 55.02 55.65 54.53 54.56 38246200    54.56
##  3     FB 2014-01-06 54.42 57.26 54.05 57.20 68852600    57.20
##  4     FB 2014-01-07 57.70 58.55 57.22 57.92 77207400    57.92
##  5     FB 2014-01-08 57.60 58.41 57.23 58.23 56682400    58.23
##  6     FB 2014-01-09 58.65 58.96 56.65 57.22 92253300    57.22
##  7     FB 2014-01-10 57.13 58.30 57.06 57.94 42449500    57.94
##  8     FB 2014-01-13 57.91 58.25 55.38 55.91 63010900    55.91
##  9     FB 2014-01-14 56.46 57.78 56.10 57.74 37503600    57.74
## 10     FB 2014-01-15 57.98 58.57 57.27 57.60 33663400    57.60
## # ... with 3,014 more rows

我们设计了一个函数来按股票代码分块绘图,可以在本文中重复使用。没有必要深究这些代码,只要认识到我们正在创建一个 ggplot2 对象,它通过指定数据框、x、y 和 group(如果存在)等要素来创建根据“symbol”分块的信息图。

# Setup plotting function that can be reused later
ggplot_facet_by_symbol <- function(data,
                                   mapping)
{
    if (is.null(mapping$group))
    {
        # No groups
        g <- data %>%
            ggplot(
                mapping = mapping)  
            labs(x = quo_name(mapping$x),
                 y = quo_name(mapping$y))
    }
    else
    {
        # Deal with groups
        g <- data %>%
            ggplot(
                mapping = mapping)   
            labs(x = quo_name(mapping$x),
                 y = quo_name(mapping$y),
                 group = quo_name(mapping$group))
    }

    # Add faceting and theme
    g <- g  
        geom_line()  
        facet_wrap(
            ~ symbol, ncol = 2, scales = "free_y")  
        scale_color_tq()  
        theme_tq()

    return(g)
}

我们可以使用绘图函数 ggplot_facet_by_symbol 快速可视化我们的数据。让我们看一下“除权调整的”股票价格。

# Plot adjusted vs date
FANG_tbl_d %>%
    ggplot_facet_by_symbol(
        mapping = aes(
            x = date, y = adjusted, color = symbol))  
    labs(
        title = "FANG Stocks: Adjusted Prices 2014 through 2016")

澳门新浦京娱乐场网站 5

上图所显示就是我们要处理的数据,下面让我们进入 tibbletime 的教程。

数据

我们将使用 tidyquant 中的 tq_get() 函数从 FRED 获取数据——啤酒、葡萄酒和蒸馏酒销售数据。

# Beer, Wine, Distilled Alcoholic Beverages, in Millions USDbeer_sales_tbl <- tq_get(    "S4248SM144NCEN",    get = "economic.data",    from = "2010-01-01",    to = "2016-12-31")beer_sales_tbl

## # A tibble: 84 x 2##          date price##        <date> <int>##  1 2010-01-01  6558##  2 2010-02-01  7481##  3 2010-03-01  9475##  4 2010-04-01  9424##  5 2010-05-01  9351##  6 2010-06-01 10552##  7 2010-07-01  9077##  8 2010-08-01  9273##  9 2010-09-01  9420## 10 2010-10-01  9413## # ... with 74 more rows

可视化数据是一个好东西,这有助于帮助我们了解正在使用的是什么数据。可视化对于时间序列分析和预测尤为重要。我们将使用 tidyquant 画图工具:主要是用 geom_ma(ma_fun = SMA,n = 12) 来添加一个周期为 12 的简单移动平均线来了解趋势。我们还可以看到似乎同时存在着趋势性(移动平均线以近似线性的模式增长)和季节性(波峰和波谷倾向于在特定月份发生)。

# Plot Beer Salesbeer_sales_tbl %>%    ggplot(aes(date, price))      geom_line(col = palette_light      geom_point(col = palette_light      geom_ma(ma_fun = SMA, n = 12, size = 1)      theme_tq()      scale_x_date(date_breaks = "1 year", date_labels = "%Y")      labs(title = "Beer Sales: 2007 through 2016")

澳门新浦京娱乐场网站 6

现在你对我们要分析的时间序列有了直观的感受,那么让我们继续!

数据

我们使用 timetk 教程中数据——啤酒、红酒和蒸馏酒销售数据,用 tidyquant 中的 tq_get() 函数从 FRED 获取。

# Beer, Wine, Distilled Alcoholic Beverages, in Millions USDbeer_sales_tbl <- tq_get(    "S4248SM144NCEN",    get = "economic.data",    from = "2010-01-01",    to = "2016-12-31")beer_sales_tbl

## # A tibble: 84 x 2##          date price##        <date> <int>##  1 2010-01-01  6558##  2 2010-02-01  7481##  3 2010-03-01  9475##  4 2010-04-01  9424##  5 2010-05-01  9351##  6 2010-06-01 10552##  7 2010-07-01  9077##  8 2010-08-01  9273##  9 2010-09-01  9420## 10 2010-10-01  9413## # ... with 74 more rows

可视化数据是一个好东西,这有助于帮助我们了解正在使用的是什么数据。可视化对于时间序列分析和预测尤为重要。我们将使用 tidyquant 画图工具:主要是用 geom_ma(ma_fun = SMA,n = 12) 来添加一个周期为 12 的简单移动平均线来了解趋势。我们还可以看到似乎同时存在着趋势性(移动平均线以近似线性的模式增长)和季节性(波峰和波谷倾向于在特定月份发生)。

# Plot Beer Salesbeer_sales_tbl %>%    ggplot(aes(date, price))      geom_line(col = palette_light      geom_point(col = palette_light      geom_ma(ma_fun = SMA, n = 12, size = 1)      theme_tq()      scale_x_date(date_breaks = "1 year", date_labels = "%Y")      labs(title = "Beer Sales: 2007 through 2016")

澳门新浦京娱乐场网站 7

现在你对我们要分析的时间序列有了直观的感受,那么让我们继续!

从 Yahoo! Finance 获得股票数据

将一列股票代码传入 tq_get(),同时设置 get = "stock.prices"。可以添加 fromto 参数设置数据的起始和结束日期。

# Get Stock Prices from Yahoo! Finance# Create a vector of stock symbolsFANG_symbols <- c("FB", "AMZN", "NFLX", "GOOG")# Pass symbols to tq_get to get daily pricesFANG_data_d <- FANG_symbols %>%    tq_get(        get = "stock.prices",        from = "2014-01-01", to = "2016-12-31")# Show the resultFANG_data_d

## # A tibble: 3,024 x 8##    symbol       date  open  high   low close   volume adjusted##     <chr>     <date> <dbl> <dbl> <dbl> <dbl>    <dbl>    <dbl>##  1     FB 2014-01-02 54.83 55.22 54.19 54.71 43195500    54.71##  2     FB 2014-01-03 55.02 55.65 54.53 54.56 38246200    54.56##  3     FB 2014-01-06 54.42 57.26 54.05 57.20 68852600    57.20##  4     FB 2014-01-07 57.70 58.55 57.22 57.92 77207400    57.92##  5     FB 2014-01-08 57.60 58.41 57.23 58.23 56682400    58.23##  6     FB 2014-01-09 58.65 58.96 56.65 57.22 92253300    57.22##  7     FB 2014-01-10 57.13 58.30 57.06 57.94 42449500    57.94##  8     FB 2014-01-13 57.91 58.25 55.38 55.91 63010900    55.91##  9     FB 2014-01-14 56.46 57.78 56.10 57.74 37503600    57.74## 10     FB 2014-01-15 57.98 58.57 57.27 57.60 33663400    57.60## # ... with 3,014 more rows

可以使用 ggplot2 画出上述结果。使用 tidyquant 提供的主题(调用 theme_tq()scale_color_tq实现金融、商务风格的可视化效果。

# Plot dataFANG_data_d %>%    ggplot(aes(x = date, y = adjusted, color = symbol))       geom_line()      facet_wrap(~ symbol, ncol = 2, scales = "free_y")      theme_tq()      scale_color_tq()      labs(title = "Visualize Financial Data")

澳门新浦京娱乐场网站 8

提示:pandas也支持基于timedeltas的指数,它可以有效代表实验或经过的时间。这本书不涉及timedelta指数,但你可以学习pandas的文档(http://pandas.pydata.org/)。

教程:tibbletime

本教程将介绍下列函数的用法:

  • filter_time:对时间索引的过滤
  • as_period:改变数据的周期
  • rollify:将任意函数转换成为滚动函数

timetk 教程:

教程分为两部分。首先,我们将遵循时间序列机器学习的工作流程。其次,我们将介绍数据转换工具。

Demo:forecast sweep 的简化预测工作流

我们将联合使用 forecastsweep 来简化预测分析。

关键想法:使用 forecast 包做预测涉及到 ts 对象,用起来并不简洁。对于 stats 包来说有 broom 来简化使用;forecast 包就用 sweep

目标:我们将用 ARIMA 模型预测未来 12 个月的数据。

从 FRED 获得经济数据

下面的例子来自房地美副首席经济学家 Leonard Kieffer 近期的文章——《A (TIDYQUANT)UM OF SOLACE》。我们将使用 tq_get() 并设置参数 get = "economic.data" 来从 FRED 经济数据库获取数据。

将一列 FRED 代码传递到 tq_get()

# Economic Data from the FRED# Create a vector of FRED symbolsFRED_symbols <- c('ETOTALUSQ176N',    # All housing units                  'EVACANTUSQ176N',   # Vacant                  'EYRVACUSQ176N',    # Year-round vacant                  'ERENTUSQ176N')     # Vacant for rent# Pass symbols to tq_get to get economic dataFRED_data_m <- FRED_symbols %>%    tq_get(get="economic.data", from = "2001-04-01")# Show resultsFRED_data_m

## # A tibble: 260 x 3##           symbol       date  price##            <chr>     <date>  <int>##  1 ETOTALUSQ176N 2001-04-01 117786##  2 ETOTALUSQ176N 2001-07-01 118216##  3 ETOTALUSQ176N 2001-10-01 118635##  4 ETOTALUSQ176N 2002-01-01 119061##  5 ETOTALUSQ176N 2002-04-01 119483##  6 ETOTALUSQ176N 2002-07-01 119909##  7 ETOTALUSQ176N 2002-10-01 120350##  8 ETOTALUSQ176N 2003-01-01 120792##  9 ETOTALUSQ176N 2003-04-01 121233## 10 ETOTALUSQ176N 2003-07-01 121682## # ... with 250 more rows

和金融数据一样,使用 ggplot2 画图,使用 tidyquant 提供的主题(调用 theme_tq()scale_color_tq实现金融、商务风格的可视化效果。

# Plot dataFRED_data_m %>%    ggplot(aes(x = date, y = price, color = symbol))       geom_line()      facet_wrap(~ symbol, ncol = 2, scales = "free_y")      theme_tq()      scale_color_tq()      labs(title = "Visualize Economic Data")

澳门新浦京娱乐场网站 9

pandas提供了许多内置的时间序列处理工具和数据算法。因此,你可以高效处理非常大的时间序列,轻松地进行切片/切块、聚合、对定期/不定期的时间序列进行重采样等。有些工具特别适合金融和经济应用,你当然也可以用它们来分析服务器日志数据。

初始化一个 tbl_time 对象

在我们使用这些新函数之前,我们需要创建一个 tbl_time 对象。新类的操作几乎与普通的 tibble 对象相同。然而,它会在背后自动跟踪时间信息。

使用 as_tbl_time() 函数初始化对象。指定 index = date,这告诉 tbl_time 对象要跟踪哪个索引。

# Convert to tbl_time
FANG_tbl_time_d <- FANG_tbl_d %>%
    as_tbl_time(index = date) 

我们可以打印 tbl_time 对象。看起来几乎与分组的 tibble 相同。请注意,“Index: date”通知我们“time tibble”已正确初始化。

# Show the tbl_time object we created
FANG_tbl_time_d

## # A time tibble: 3,024 x 8
## # Index:  date
## # Groups: symbol [4]
##    symbol       date  open  high   low close   volume adjusted
##     <chr>     <date> <dbl> <dbl> <dbl> <dbl>    <dbl>    <dbl>
##  1     FB 2014-01-02 54.83 55.22 54.19 54.71 43195500    54.71
##  2     FB 2014-01-03 55.02 55.65 54.53 54.56 38246200    54.56
##  3     FB 2014-01-06 54.42 57.26 54.05 57.20 68852600    57.20
##  4     FB 2014-01-07 57.70 58.55 57.22 57.92 77207400    57.92
##  5     FB 2014-01-08 57.60 58.41 57.23 58.23 56682400    58.23
##  6     FB 2014-01-09 58.65 58.96 56.65 57.22 92253300    57.22
##  7     FB 2014-01-10 57.13 58.30 57.06 57.94 42449500    57.94
##  8     FB 2014-01-13 57.91 58.25 55.38 55.91 63010900    55.91
##  9     FB 2014-01-14 56.46 57.78 56.10 57.74 37503600    57.74
## 10     FB 2014-01-15 57.98 58.57 57.27 57.60 33663400    57.60
## # ... with 3,014 more rows

我们可以使用绘图函数 ggplot_facet_by_symbol() 绘制它,我们看到 tbl_time 对象与 tbl 对象的反应相同。

# Plot the tbl_time object
FANG_tbl_time_d %>%
    ggplot_facet_by_symbol(
        mapping = aes(
            x = date, y = adjusted, color = symbol))  
    labs(
        title = "Working with tbltime: Reacts same as tbl class")

澳门新浦京娱乐场网站 10

PART 1:时间序列机器学习

时间序列机器学习是预测时间序列数据的一种很好的方法,但在我们开始之前,这里有几个注意点:

  • 关键洞察力:将时间序列签名(时间戳信息按列扩展到特征集)用于执行机器学习。
  • 目标:我们将使用时间序列签名预测未来 12 个月的时间序列数据。

我们将遵循可用于执行时间序列机器学习的工作流程。你将看到几个 timetk 函数如何帮助完成此过程。我们将使用简单的 lm() 线性回归进行机器学习,你将看到使用时间序列签名会使机器学习更强大和准确。此外,你还应该考虑使用其他更强大的机器学习算法,例如 xgboostglmnet等。

STEP 1:创建 ts 对象

使用 timetk::tk_ts()tbl 转换成 ts,从之前的教程可以了解到这个函数有两点好处:

  1. 这是一个统一的方法,实现与 ts 对象的相互转换。
  2. 得到的 ts 对象包含 timetk_idx 属性,是一个基于初始时间信息的索引。

下面开始转换,注意 ts 对象是规则时间序列,所以要设置 startfreq

# Convert from tbl to tsbeer_sales_ts <- tk_ts(    beer_sales_tbl,    start = 2010,    freq = 12)beer_sales_ts

##        Jan   Feb   Mar   Apr   May   Jun   Jul   Aug   Sep   Oct## 2010  6558  7481  9475  9424  9351 10552  9077  9273  9420  9413## 2011  6901  8014  9833  9281  9967 11344  9106 10468 10085  9612## 2012  7486  8641  9709  9423 11342 11274  9845 11163  9532 10754## 2013  8395  8888 10109 10493 12217 11385 11186 11462 10494 11541## 2014  8559  9061 10058 10979 11794 11906 10966 10981 10827 11815## 2015  8398  9061 10720 11105 11505 12903 11866 11223 12023 11986## 2016  8540 10158 11879 11155 11916 13291 10540 12212 11786 11424##        Nov   Dec## 2010  9866 11455## 2011 10328 11483## 2012 10953 11922## 2013 11139 12709## 2014 10466 13303## 2015 11510 14190## 2016 12482 13832

检查 ts 对象具有 timetk_idx 属性。

# Check that ts-object has a timetk indexhas_timetk_idx(beer_sales_ts)

## [1] TRUE

OK,这对后面要用的 sw_sweep() 很重要。下面我们就要建立 ARIMA 模型了。

使用 tq_transmutetq_mutate 转换数据

函数 tq_transmute()tq_mutate() 可以使 xtszooquantmod 中的函数调用更“tidy”。这里直接介绍使用,“可用函数”一节罗列了已经整合进 tidyquant 的若干其他函数。

11.1 日期和时间数据类型及工具

Python标准库包含用于日期(date)和时间(time)数据的数据类型,而且还有日历方面的功能。我们主要会用到datetime、time以及calendar模块。datetime.datetime(也可以简写为datetime)是用得最多的数据类型:

In [10]: from datetime import datetime

In [11]: now = datetime.now()

In [12]: now
Out[12]: datetime.datetime(2017, 9, 25, 14, 5, 52, 72973)

In [13]: now.year, now.month, now.day
Out[13]: (2017, 9, 25)

datetime以毫秒形式存储日期和时间。timedelta表示两个datetime对象之间的时间差:

In [14]: delta = datetime(2011, 1, 7) - datetime(2008, 6, 24, 8, 15)

In [15]: delta
Out[15]: datetime.timedelta(926, 56700)

In [16]: delta.days
Out[16]: 926

In [17]: delta.seconds
Out[17]: 56700

可以给datetime对象加上(或减去)一个或多个timedelta,这样会产生一个新对象:

In [18]: from datetime import timedelta

In [19]: start = datetime(2011, 1, 7)

In [20]: start   timedelta(12)
Out[20]: datetime.datetime(2011, 1, 19, 0, 0)

In [21]: start - 2 * timedelta(12)
Out[21]: datetime.datetime(2010, 12, 14, 0, 0)

datetime模块中的数据类型参见表10-1。虽然本章主要讲的是pandas数据类型和高级时间序列处理,但你肯定会在Python的其他地方遇到有关datetime的数据类型。

表11-1 datetime模块中的数据类型

澳门新浦京娱乐场网站 11

tzinfo 存储时区信息的基本类型

时间序列函数

让我们看看可以用新的 tbl_time 对象做些什么。

STEP 0:回顾数据

看看我们的起点,先打印出 beer_sales_tbl

# Starting pointbeer_sales_tbl

## # A tibble: 84 x 2##          date price##        <date> <int>##  1 2010-01-01  6558##  2 2010-02-01  7481##  3 2010-03-01  9475##  4 2010-04-01  9424##  5 2010-05-01  9351##  6 2010-06-01 10552##  7 2010-07-01  9077##  8 2010-08-01  9273##  9 2010-09-01  9420## 10 2010-10-01  9413## # ... with 74 more rows

我们可以使用 tk_index() 来提取索引;使用 tk_get_timeseries_summary() 来检索索引的摘要信息,进而快速了解时间序列。我们使用 glimpse() 输出一个更漂亮的格式。

beer_sales_tbl %>%    tk_index() %>%    tk_get_timeseries_summary() %>%    glimpse()

## Observations: 1## Variables: 12## $ n.obs        <int> 84## $ start        <date> 2010-01-01## $ end          <date> 2016-12-01## $ units        <chr> "days"## $ scale        <chr> "month"## $ tzone        <chr> "UTC"## $ diff.minimum <dbl> 2419200## $ diff.q1      <dbl> 2592000## $ diff.median  <dbl> 2678400## $ diff.mean    <dbl> 2629475## $ diff.q3      <dbl> 2678400## $ diff.maximum <dbl> 2678400

你可以看到一些重要的特征,例如起始、结束、单位等等。还有时间差的分位数(相邻两个观察之间差距的秒数),这对评估规律性的程度很有用。由于时间尺度是月度的,因此每个月之间差距的秒数并不规则。

STEP 2A:ARIMA 模型

我们使用 forecast 包里的 auto.arima() 函数为时间序列建模。

# Model using auto.arimafit_arima <- auto.arima(beer_sales_ts)fit_arima

## Series: beer_sales_ts ## ARIMA[12] with drift ## ## Coefficients:##           ar1     ar2     ar3     sar1    drift##       -0.2498  0.1079  0.6210  -0.2817  32.1157## s.e.   0.0933  0.0982  0.0925   0.1333   5.8882## ## sigma^2 estimated as 175282:  log likelihood=-535.49## AIC=1082.97   AICc=1084.27   BIC=1096.63

tq_transmute

tq_transmute()tq_mutate() 之间的区别在于 tq_transmute() 将返回一个新的数据框对象,而 tq_mutate() 则在原有数据框的基础上横向添加数据。当数据因为改变周期而改变行数时,tq_transmute() 特别有用。

字符串和datetime的相互转换

利用str或strftime方法(传入一个格式化字符串),datetime对象和pandas的Timestamp对象(稍后就会介绍)可以被格式化为字符串:

In [22]: stamp = datetime(2011, 1, 3)

In [23]: str(stamp)
Out[23]: '2011-01-03 00:00:00'

In [24]: stamp.strftime('%Y-%m-%d')
Out[24]: '2011-01-03'

表11-2列出了全部的格式化编码。

表11-2 datetime格式定义(兼容ISO C89)

澳门新浦京娱乐场网站 12

澳门新浦京娱乐场网站 13

datetime.strptime可以用这些格式化编码将字符串转换为日期:

In [25]: value = '2011-01-03'

In [26]: datetime.strptime(value, '%Y-%m-%d')
Out[26]: datetime.datetime(2011, 1, 3, 0, 0)

In [27]: datestrs = ['7/6/2011', '8/6/2011']

In [28]: [datetime.strptime(x, '%m/%d/%Y') for x in datestrs]
Out[28]: 
[datetime.datetime(2011, 7, 6, 0, 0),
 datetime.datetime(2011, 8, 6, 0, 0)]

datetime.strptime是通过已知格式进行日期解析的最佳方式。但是每次都要编写格式定义是很麻烦的事情,尤其是对于一些常见的日期格式。这种情况下,你可以用dateutil这个第三方包中的parser.parse方法(pandas中已经自动安装好了):

In [29]: from dateutil.parser import parse

In [30]: parse('2011-01-03')
Out[30]: datetime.datetime(2011, 1, 3, 0, 0)

dateutil可以解析几乎所有人类能够理解的日期表示形式:

In [31]: parse('Jan 31, 1997 10:45 PM')
Out[31]: datetime.datetime(1997, 1, 31, 22, 45)

在国际通用的格式中,日出现在月的前面很普遍,传入dayfirst=True即可解决这个问题:

In [32]: parse('6/12/2011', dayfirst=True)
Out[32]: datetime.datetime(2011, 12, 6, 0, 0)

pandas通常是用于处理成组日期的,不管这些日期是DataFrame的轴索引还是列。to_datetime方法可以解析多种不同的日期表示形式。对标准日期格式(如ISO8601)的解析非常快:

In [33]: datestrs = ['2011-07-06 12:00:00', '2011-08-06 00:00:00']

In [34]: pd.to_datetime(datestrs)
Out[34]: DatetimeIndex(['2011-07-06 12:00:00', '2011-08-06 00:00:00'], dtype='dat
etime64[ns]', freq=None)

它还可以处理缺失值(None、空字符串等):

In [35]: idx = pd.to_datetime(datestrs   [None])

In [36]: idx
Out[36]: DatetimeIndex(['2011-07-06 12:00:00', '2011-08-06 00:00:00', 'NaT'], dty
pe='datetime64[ns]', freq=None)

In [37]: idx[2]
Out[37]: NaT

In [38]: pd.isnull(idx)
Out[38]: array([False, False,  True], dtype=bool)

NaT(Not a Time)是pandas中时间戳数据的null值。

注意:dateutil.parser是一个实用但不完美的工具。比如说,它会把一些原本不是日期的字符串认作是日期(比如"42"会被解析为2042年的今天)。

datetime对象还有一些特定于当前环境(位于不同国家或使用不同语言的系统)的格式化选项。例如,德语或法语系统所用的月份简写就与英语系统所用的不同。表11-3进行了总结。

表11-3 特定于当前环境的日期格式

澳门新浦京娱乐场网站 14

filter_time

filter_time() 函数根据按日期简便快捷地过滤 tbl_time 对象,它使用一个函数格式(例如 'date_operator_start' ~ 'date_operator_end')。我们使用标准日期格式 YYYY-MM-DD HH:MM:SS 指定日期运算符,但也有强大的简化标记来更有效地指定日期子集。

假设我们想要过滤出 2014-06-012014-06-15 之间的所有观察结果。我们可以使用函数标记 filter_time('2014-06-01' ~ '2014-06-15') 来完成。

# filter_time by day
FANG_tbl_time_d %>%
    filter_time('2014-06-01' ~ '2014-06-15') %>%
    # Plotting
    ggplot_facet_by_symbol(
        mapping = aes(
            x = date, y = adjusted, color = symbol))  
    geom_point()  
    labs(
        title = "Time Filter: Use functional notation to quickly subset by time",
        subtitle = "2014-06-01 ~ 2014-06-15")

澳门新浦京娱乐场网站 15

我们可以按月完成同样的工作。假设我们只想在 2014 年 3 月进行观察。使用简化函数标记 ~ '2014-03'

# filter_time by month
FANG_tbl_time_d %>%
    filter_time(~ '2014-03') %>%
    # Plotting
    ggplot_facet_by_symbol(
        mapping = aes(
            x = date, y = adjusted, color = symbol))  
    geom_point()  
    labs(
        title = "Time Filter: Use shorthand for even easier subsetting",
        subtitle = "~ 2014-03")

澳门新浦京娱乐场网站 16

tbl_time 对象也响应括号符号运算符——[。在这里,我们提取 2014 年所有日期的数据。

# time filter bracket [] notation
FANG_tbl_time_d[~ '2014'] %>%
    # Plotting
    ggplot_facet_by_symbol(
        mapping = aes(
            x = date, y = adjusted, color = symbol))  
    labs(
        title = "Time Filter: Bracket Notation Works Too",
        subtitle = "FANG_tbl_time_d[~ 2014]")

澳门新浦京娱乐场网站 17

filter_time() 有许多功能和简化标记,感兴趣的读者可以查看 filter_time vignette 和 filter_time function documentation。

STEP 1:扩充时间序列签名

tk_augment_timeseries_signature() 函数将时间戳信息逐列扩展到机器学习特征集中,并将时间序列信息列添加到初始数据表。

# Augment (adds data frame columns)beer_sales_tbl_aug <- beer_sales_tbl %>%    tk_augment_timeseries_signature()beer_sales_tbl_aug

## # A tibble: 84 x 30##          date price  index.num    diff  year year.iso  half quarter##        <date> <int>      <int>   <int> <int>    <int> <int>   <int>##  1 2010-01-01  6558 1262304000      NA  2010     2009     1       1##  2 2010-02-01  7481 1264982400 2678400  2010     2010     1       1##  3 2010-03-01  9475 1267401600 2419200  2010     2010     1       1##  4 2010-04-01  9424 1270080000 2678400  2010     2010     1       2##  5 2010-05-01  9351 1272672000 2592000  2010     2010     1       2##  6 2010-06-01 10552 1275350400 2678400  2010     2010     1       2##  7 2010-07-01  9077 1277942400 2592000  2010     2010     2       3##  8 2010-08-01  9273 1280620800 2678400  2010     2010     2       3##  9 2010-09-01  9420 1283299200 2678400  2010     2010     2       3## 10 2010-10-01  9413 1285891200 2592000  2010     2010     2       4## # ... with 74 more rows, and 22 more variables: month <int>,## #   month.xts <int>, month.lbl <ord>, day <int>, hour <int>,## #   minute <int>, second <int>, hour12 <int>, am.pm <int>,## #   wday <int>, wday.xts <int>, wday.lbl <ord>, mday <int>,## #   qday <int>, yday <int>, mweek <int>, week <int>, week.iso <int>,## #   week2 <int>, week3 <int>, week4 <int>, mday7 <int>

STEP 2B:简化模型

就像 broom 简化 stats 包的使用一样,我么可以使用 sweep 的函数简化 ARIMA 模型。下面介绍三个函数:

  • sw_tidy():用于检索模型参数
  • sw_glance():用于检索模型描述和训练集的精确度度量
  • sw_augment():用于获得模型残差

改变周期与 tq_transmute

下面的例子将改变数据的周期,从每日数据变为月度数据。这时,你需要使用 tq_transmute() 来完成这一操作,因为数据的行数改变了。

# Change periodicity from daily to monthly using to.period from xtsFANG_data_m <- FANG_data_d %>%    group_by %>%    tq_transmute(        select      = adjusted,        mutate_fun  = to.period,        period      = "months")FANG_data_m

## # A tibble: 144 x 3## # Groups:   symbol [4]##    symbol       date adjusted##     <chr>     <date>    <dbl>##  1     FB 2014-01-31    62.57##  2     FB 2014-02-28    68.46##  3     FB 2014-03-31    60.24##  4     FB 2014-04-30    59.78##  5     FB 2014-05-30    63.30##  6     FB 2014-06-30    67.29##  7     FB 2014-07-31    72.65##  8     FB 2014-08-29    74.82##  9     FB 2014-09-30    79.04## 10     FB 2014-10-31    74.99## # ... with 134 more rows

改变数据周期可以缩减数据量。一些注意事项项:

  • theme_tq()scale_color_tq() 用来绘制商务风格的图。
  • 如果要改变数据的周期或者进行其他基于时间的操作,请留意后续关于 tibbletime 包的教程,tibbletime 以另外一种标准处理基于时间的操作。

11.2 时间序列基础

pandas最基本的时间序列类型就是以时间戳(通常以Python字符串或datatime对象表示)为索引的Series:

In [39]: from datetime import datetime

In [40]: dates = [datetime(2011, 1, 2), datetime(2011, 1, 5),
   ....:          datetime(2011, 1, 7), datetime(2011, 1, 8),
   ....:          datetime(2011, 1, 10), datetime(2011, 1, 12)]

In [41]: ts = pd.Series(np.random.randn(6), index=dates)

In [42]: ts
Out[42]: 
2011-01-02   -0.204708
2011-01-05    0.478943
2011-01-07   -0.519439
2011-01-08   -0.555730
2011-01-10    1.965781
2011-01-12    1.393406
dtype: float64

这些datetime对象实际上是被放在一个DatetimeIndex中的:

In [43]: ts.index
Out[43]: 
DatetimeIndex(['2011-01-02', '2011-01-05', '2011-01-07', '2011-01-08',
               '2011-01-10', '2011-01-12'],
              dtype='datetime64[ns]', freq=None)

跟其他Series一样,不同索引的时间序列之间的算术运算会自动按日期对齐:

In [44]: ts   ts[::2]
Out[44]: 
2011-01-02   -0.409415
2011-01-05         NaN
2011-01-07   -1.038877
2011-01-08         NaN
2011-01-10    3.931561
2011-01-12         NaN
dtype: float64

ts[::2] 是每隔两个取一个。

pandas用NumPy的datetime64数据类型以纳秒形式存储时间戳:

In [45]: ts.index.dtype
Out[45]: dtype('<M8[ns]')

DatetimeIndex中的各个标量值是pandas的Timestamp对象:

In [46]: stamp = ts.index[0]

In [47]: stamp
Out[47]: Timestamp('2011-01-02 00:00:00')

只要有需要,TimeStamp可以随时自动转换为datetime对象。此外,它还可以存储频率信息(如果有的话),且知道如何执行时区转换以及其他操作。稍后将对此进行详细讲解。

as_period

函数 as_period() 可以改变 tbl_time 对象的周期。与传统方法相比,使用此方法有两个优点:

  1. 函数标记非常灵活:yearly == y == 1 y
  2. 函数标记提供了无数周期转换的可能,例如:
    • 15 d:以 15 天为一周期
    • 2 m:以 2 月为一周期
    • 4 m:以 4 月为一周期
    • 6 m:以半年为一周期

首先,让我们做一个简单的月度周期性变化。

# Convert from daily to monthly periodicity
FANG_tbl_time_d %>%
    as_period(period = "month") %>%
    # Plotting
    ggplot_facet_by_symbol(
        mapping = aes(
            x = date, y = adjusted, color = symbol))  
    labs(
        title = "Periodicity Change from Daily to Monthly")  
    geom_point()

澳门新浦京娱乐场网站 18

让我们提升一个档次。那么每两个月一次呢? 只需使用函数标记 2 m 即可。

# Convert from daily to bi-monthly periodicity
FANG_tbl_time_d %>%
    as_period(period = '2 m') %>%
    # Plotting
    ggplot_facet_by_symbol(
        mapping = aes(
            x = date, y = adjusted, color = symbol))  
    labs(
        title = "Periodicity Change to Daily to Bi-Monthly",
        subtitle = "2~m")  
    geom_point()

澳门新浦京娱乐场网站 19

让我们继续。那么每半年一次呢? 只需使用 6 m 即可。

# Convert from daily to bi-annually periodicity
FANG_tbl_time_d %>%
    as_period(period = '6 m') %>%
    # Plotting
    ggplot_facet_by_symbol(
        mapping = aes(
            x = date, y = adjusted, color = symbol))  
    labs(
        title = "Periodicity Change to Daily to Bi-Annually",
        subtitle = "6~m")  
    geom_point()

澳门新浦京娱乐场网站 20

函数标记几乎提供了无限可能,感兴趣的话可以查看 vignette on periodicity change with tibbletime。

STEP 2:模型

任何回归模型都可以应用于数据,我们在这里使用 lm()。 请注意,我们删除了 datediff 列。大多数算法无法使用日期数据,而 diff 列对机器学习没有什么用处(它对于查找数据中的时间间隔更有用)。

# linear regression model used, but can use any modelfit_lm <- lm(    price ~ .,    data = select(        beer_sales_tbl_aug,        -c(date, diff)))summary

## ## Call:## lm(formula = price ~ ., data = select(beer_sales_tbl_aug, -c(date, ##     diff)))## ## Residuals:##    Min     1Q Median     3Q    Max ## -447.3 -145.4  -18.2  169.8  421.4 ## ## Coefficients: (16 not defined because of singularities)##                Estimate Std. Error t value Pr    ## (Intercept)   3.660e 08  1.245e 08   2.940 0.004738 ** ## index.num     5.900e-03  2.003e-03   2.946 0.004661 ** ## year         -1.974e 05  6.221e 04  -3.173 0.002434 ** ## year.iso      1.159e 04  6.546e 03   1.770 0.082006 .  ## half         -2.132e 03  6.107e 02  -3.491 0.000935 ***## quarter      -1.239e 04  2.190e 04  -0.566 0.573919    ## month        -3.910e 03  7.355e 03  -0.532 0.597058    ## month.xts            NA         NA      NA       NA    ## month.lbl.L          NA         NA      NA       NA    ## month.lbl.Q  -1.643e 03  2.069e 02  -7.942 8.59e-11 ***## month.lbl.C   8.368e 02  5.139e 02   1.628 0.108949    ## month.lbl^4   6.452e 02  1.344e 02   4.801 1.18e-05 ***## month.lbl^5   7.563e 02  4.241e 02   1.783 0.079852 .  ## month.lbl^6   3.206e 02  1.609e 02   1.992 0.051135 .  ## month.lbl^7  -3.537e 02  1.867e 02  -1.894 0.063263 .  ## month.lbl^8   3.687e 02  3.217e 02   1.146 0.256510    ## month.lbl^9          NA         NA      NA       NA    ## month.lbl^10  6.339e 02  2.240e 02   2.830 0.006414 ** ## month.lbl^11         NA         NA      NA       NA    ## day                  NA         NA      NA       NA    ## hour                 NA         NA      NA       NA    ## minute               NA         NA      NA       NA    ## second               NA         NA      NA       NA    ## hour12               NA         NA      NA       NA    ## am.pm                NA         NA      NA       NA    ## wday         -8.264e 01  1.898e 01  -4.353 5.63e-05 ***## wday.xts             NA         NA      NA       NA    ## wday.lbl.L           NA         NA      NA       NA    ## wday.lbl.Q   -7.109e 02  1.093e 02  -6.503 2.13e-08 ***## wday.lbl.C    2.355e 02  1.336e 02   1.763 0.083273 .  ## wday.lbl^4    8.033e 01  1.133e 02   0.709 0.481281    ## wday.lbl^5    6.480e 01  8.029e 01   0.807 0.422951    ## wday.lbl^6    2.276e 01  8.200e 01   0.278 0.782319    ## mday                 NA         NA      NA       NA    ## qday         -1.362e 02  2.418e 02  -0.563 0.575326    ## yday         -2.356e 02  1.416e 02  -1.664 0.101627    ## mweek        -1.670e 02  1.477e 02  -1.131 0.262923    ## week         -1.764e 02  1.890e 02  -0.933 0.354618    ## week.iso      2.315e 02  1.256e 02   1.842 0.070613 .  ## week2        -1.235e 02  1.547e 02  -0.798 0.428283    ## week3                NA         NA      NA       NA    ## week4                NA         NA      NA       NA    ## mday7                NA         NA      NA       NA    ## ---## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1## ## Residual standard error: 260.4 on 57 degrees of freedom## Multiple R-squared:  0.9798, Adjusted R-squared:  0.9706 ## F-statistic: 106.4 on 26 and 57 DF,  p-value: < 2.2e-16

sw_tidy

sw_tidy() 函数以 tibble 对象的形式返回模型参数。

# sw_tidy - Get model coefficientssw_tidy(fit_arima)

## # A tibble: 5 x 2##    term   estimate##   <chr>      <dbl>## 1   ar1 -0.2497937## 2   ar2  0.1079269## 3   ar3  0.6210345## 4  sar1 -0.2816877## 5 drift 32.1157478

周期改变前,数据太多

# Daily dataFANG_data_d %>%    ggplot(aes(date, adjusted, color = symbol))      geom_point()      geom_line()      facet_wrap(~ symbol, ncol = 2, scales = "free_y")      scale_color_tq()      theme_tq()      labs(title = "Before transformation: Too Much Data")

澳门新浦京娱乐场网站 21

索引、选取、子集构造

当你根据标签索引选取数据时,时间序列和其它的pandas.Series很像:

In [48]: stamp = ts.index[2]

In [49]: ts[stamp]
Out[49]: -0.51943871505673811

还有一种更为方便的用法:传入一个可以被解释为日期的字符串:

In [50]: ts['1/10/2011']
Out[50]: 1.9657805725027142

In [51]: ts['20110110']
Out[51]: 1.9657805725027142

对于较长的时间序列,只需传入“年”或“年月”即可轻松选取数据的切片:

In [52]: longer_ts = pd.Series(np.random.randn(1000),
   ....:                       index=pd.date_range('1/1/2000', periods=1000))

In [53]: longer_ts
Out[53]: 
2000-01-01    0.092908
2000-01-02    0.281746
2000-01-03    0.769023
2000-01-04    1.246435
2000-01-05    1.007189
2000-01-06   -1.296221
2000-01-07    0.274992
2000-01-08    0.228913
2000-01-09    1.352917
2000-01-10    0.886429
                ...   
2002-09-17   -0.139298
2002-09-18   -1.159926
2002-09-19    0.618965
2002-09-20    1.373890
2002-09-21   -0.983505
2002-09-22    0.930944
2002-09-23   -0.811676
2002-09-24   -1.830156
2002-09-25   -0.138730
2002-09-26    0.334088
Freq: D, Length: 1000, dtype: float64

In [54]: longer_ts['2001']
Out[54]: 
2001-01-01    1.599534
2001-01-02    0.474071
2001-01-03    0.151326
2001-01-04   -0.542173
2001-01-05   -0.475496
2001-01-06    0.106403
2001-01-07   -1.308228
2001-01-08    2.173185
2001-01-09    0.564561
2001-01-10   -0.190481
                ...   
2001-12-22    0.000369
2001-12-23    0.900885
2001-12-24   -0.454869
2001-12-25   -0.864547
2001-12-26    1.129120
2001-12-27    0.057874
2001-12-28   -0.433739
2001-12-29    0.092698
2001-12-30   -1.397820
2001-12-31    1.457823
Freq: D, Length: 365, dtype: float64

这里,字符串“2001”被解释成年,并根据它选取时间区间。指定月也同样奏效:

In [55]: longer_ts['2001-05']
Out[55]: 
2001-05-01   -0.622547
2001-05-02    0.936289
2001-05-03    0.750018
2001-05-04   -0.056715
2001-05-05    2.300675
2001-05-06    0.569497
2001-05-07    1.489410
2001-05-08    1.264250
2001-05-09   -0.761837
2001-05-10   -0.331617
                ...   
2001-05-22    0.503699
2001-05-23   -1.387874
2001-05-24    0.204851
2001-05-25    0.603705
2001-05-26    0.545680
2001-05-27    0.235477
2001-05-28    0.111835
2001-05-29   -1.251504
2001-05-30   -2.949343
2001-05-31    0.634634
Freq: D, Length: 31, dtype: float64

datetime对象也可以进行切片:

In [56]: ts[datetime(2011, 1, 7):]
Out[56]: 
2011-01-07   -0.519439
2011-01-08   -0.555730
2011-01-10    1.965781
2011-01-12    1.393406
dtype: float64

由于大部分时间序列数据都是按照时间先后排序的,因此你也可以用不存在于该时间序列中的时间戳对其进行切片(即范围查询):

In [57]: ts
Out[57]: 
2011-01-02   -0.204708
2011-01-05    0.478943
2011-01-07   -0.519439
2011-01-08   -0.555730
2011-01-10    1.965781
2011-01-12    1.393406
dtype: float64

In [58]: ts['1/6/2011':'1/11/2011']
Out[58]: 
2011-01-07   -0.519439
2011-01-08   -0.555730
2011-01-10    1.965781
dtype: float64

跟之前一样,你可以传入字符串日期、datetime或Timestamp。注意,这样切片所产生的是源时间序列的视图,跟NumPy数组的切片运算是一样的。

这意味着,没有数据被复制,对切片进行修改会反映到原始数据上。

此外,还有一个等价的实例方法也可以截取两个日期之间TimeSeries:

In [59]: ts.truncate(after='1/9/2011')
Out[59]: 
2011-01-02   -0.204708
2011-01-05    0.478943
2011-01-07   -0.519439
2011-01-08   -0.555730
dtype: float64

面这些操作对DataFrame也有效。例如,对DataFrame的行进行索引:

In [60]: dates = pd.date_range('1/1/2000', periods=100, freq='W-WED')

In [61]: long_df = pd.DataFrame(np.random.randn(100, 4),
   ....:                        index=dates,
   ....:                        columns=['Colorado', 'Texas',
   ....:                                 'New York', 'Ohio'])

In [62]: long_df.loc['5-2001']
Out[62]: 
            Colorado     Texas  New York      Ohio
2001-05-02 -0.006045  0.490094 -0.277186 -0.707213
2001-05-09 -0.560107  2.735527  0.927335  1.513906
2001-05-16  0.538600  1.273768  0.667876 -0.969206
2001-05-23  1.676091 -0.817649  0.050188  1.951312
2001-05-30  3.260383  0.963301  1.201206 -1.852001

rollify

rollify() 函数是一个副词tidyverse 中的一种特殊类型的函数,用于修改另一个函数)。rollify() 的作用是将任何函数转换为自身的滚动版本。

# Rolling 60-day mean
roll_mean_60 <- rollify(
    mean, window = 60)

FANG_tbl_time_d %>%
    mutate(
        mean_60 = roll_mean_60(adjusted)) %>%
    select(-c(open:volume)) %>%
    # Plot
    ggplot_facet_by_symbol(
        mapping = aes(
            x = date, y = adjusted, color = symbol))  
    geom_line(
        aes(y = mean_60),
        color = palette_light()[[6]])  
    labs(
        title = "Rolling 60-Day Mean with rollify")

澳门新浦京娱乐场网站 22

我们甚至可以做出更复杂的滚动功能,例如相关性。我们在 rollify() 中使用函数形式 .f = ~fun(.x,.y,...)

# Rolling correlation
roll_corr_60 <- rollify(
    ~ cor(.x, .y, use = "pairwise.complete.obs"),
    window = 60)

FANG_tbl_time_d %>%
    mutate(
        cor_60 = roll_corr_60(
            open, close)) %>%
    select(-c(open:adjusted)) %>%
    # Plot
    ggplot_facet_by_symbol(
        mapping = aes(
            x = date, y = cor_60, color = symbol))  
    labs(
        title = "Rollify: 60-Day Rolling Correlation Between Open and Close Prices")

澳门新浦京娱乐场网站 23

我们甚至可以返回多个结果。例如,我们可以创建滚动分位数。

首先,创建一个返回分位数的函数。

# Quantile tbl function
quantile_tbl <- function(x)
{
    q <- quantile(x)

    tibble(
        quantile_name  = names(q),
        quantile_value = q)
}

# Test the function
quantile_tbl(1:100)

## # A tibble: 5 x 2
##   quantile_name quantile_value
##           <chr>          <dbl>
## 1            0%           1.00
## 2           25%          25.75
## 3           50%          50.50
## 4           75%          75.25
## 5          100%         100.00

很好,它可以工作。接下来,使用 rollify 创建滚动版本。我们设置 unlist = FALSE 来返回列表列

# Rollified quantile function
roll_quantile_60 <- rollify(
    quantile_tbl, window = 60, unlist = FALSE)

接下来,在 mutate() 中应用滚动分位数函数来获得滚动分位数。确保你已经用 select()filter()unnest() 删除了不必要的列,过滤了 NA 值,并展开列表列。现在每个日期有五个分位数值。

# Apply rolling quantile 
FANG_quantile_60 <- FANG_tbl_time_d %>%
    mutate(
        rolling_quantile = roll_quantile_60(adjusted)) %>%
    select(-c(open:adjusted)) %>%
    filter(!is.na(rolling_quantile)) %>%
    unnest()

FANG_quantile_60

## # A time tibble: 13,940 x 4
## # Index:  date
## # Groups: symbol [4]
##    symbol       date quantile_name quantile_value
##  *  <chr>     <date>         <chr>          <dbl>
##  1     FB 2014-03-28            0%        53.5300
##  2     FB 2014-03-28           25%        57.8750
##  3     FB 2014-03-28           50%        64.2100
##  4     FB 2014-03-28           75%        68.6275
##  5     FB 2014-03-28          100%        72.0300
##  6     FB 2014-03-31            0%        53.5300
##  7     FB 2014-03-31           25%        57.9350
##  8     FB 2014-03-31           50%        64.2100
##  9     FB 2014-03-31           75%        68.6275
## 10     FB 2014-03-31          100%        72.0300
## # ... with 13,930 more rows

最后,画出结果。

FANG_quantile_60 %>%
    ggplot_facet_by_symbol(
        mapping = aes(
            x = date, y = quantile_value,
            color = symbol, group = quantile_name))  
    labs(
        title = "Rollify: Create Rolling Quantiles")

澳门新浦京娱乐场网站 24

如果想继续探索 rollify 的用法,可以查看 vignette on rolling functions with rollify。

STEP 3:构建未来数据

使用 tk_index() 来扩展索引。

# Retrieves the timestamp informationbeer_sales_idx <- beer_sales_tbl %>%    tk_index()tail(beer_sales_idx)

## [1] "2016-07-01" "2016-08-01" "2016-09-01" "2016-10-01" "2016-11-01"## [6] "2016-12-01"

通过 tk_make_future_timeseries 函数从现有索引构造未来索引。函数会在内部检查索引的周期性,并返回正确的序列。我们在 whole vignette on how to make future time series 介绍了该主题更详尽的内容。

# Make future indexfuture_idx <- beer_sales_idx %>%    tk_make_future_timeseries(        n_future = 12)future_idx

##  [1] "2017-01-01" "2017-02-01" "2017-03-01" "2017-04-01" "2017-05-01"##  [6] "2017-06-01" "2017-07-01" "2017-08-01" "2017-09-01" "2017-10-01"## [11] "2017-11-01" "2017-12-01"

tk_get_timeseries_signature() 将未来索引转换成时间序列签名数据框。

new_data_tbl <- future_idx %>%    tk_get_timeseries_signature()new_data_tbl

## # A tibble: 12 x 29##         index  index.num    diff  year year.iso  half quarter month##        <date>      <int>   <int> <int>    <int> <int>   <int> <int>##  1 2017-01-01 1483228800      NA  2017     2016     1       1     1##  2 2017-02-01 1485907200 2678400  2017     2017     1       1     2##  3 2017-03-01 1488326400 2419200  2017     2017     1       1     3##  4 2017-04-01 1491004800 2678400  2017     2017     1       2     4##  5 2017-05-01 1493596800 2592000  2017     2017     1       2     5##  6 2017-06-01 1496275200 2678400  2017     2017     1       2     6##  7 2017-07-01 1498867200 2592000  2017     2017     2       3     7##  8 2017-08-01 1501545600 2678400  2017     2017     2       3     8##  9 2017-09-01 1504224000 2678400  2017     2017     2       3     9## 10 2017-10-01 1506816000 2592000  2017     2017     2       4    10## 11 2017-11-01 1509494400 2678400  2017     2017     2       4    11## 12 2017-12-01 1512086400 2592000  2017     2017     2       4    12## # ... with 21 more variables: month.xts <int>, month.lbl <ord>,## #   day <int>, hour <int>, minute <int>, second <int>, hour12 <int>,## #   am.pm <int>, wday <int>, wday.xts <int>, wday.lbl <ord>,## #   mday <int>, qday <int>, yday <int>, mweek <int>, week <int>,## #   week.iso <int>, week2 <int>, week3 <int>, week4 <int>,## #   mday7 <int>

sw_glance

sw_glance() 函数以 tibble 对象的形式返回训练集的精确度度量。可以使用 glimpse 函数美化显示结果。

# sw_glance - Get model description and training set accuracy measuressw_glance(fit_arima) %>%    glimpse()

## Observations: 1## Variables: 12## $ model.desc <chr> "ARIMA[12] with drift"## $ sigma      <dbl> 418.6665## $ logLik     <dbl> -535.4873## $ AIC        <dbl> 1082.975## $ BIC        <dbl> 1096.635## $ ME         <dbl> 1.189875## $ RMSE       <dbl> 373.9091## $ MAE        <dbl> 271.7068## $ MPE        <dbl> -0.06716239## $ MAPE       <dbl> 2.526077## $ MASE       <dbl> 0.4989005## $ ACF1       <dbl> 0.02215405

周期改变后,容易理解

tq_transmute() 转变成月度数据后容易理解多了。

# Monthly dataFANG_data_m %>%    ggplot(aes(date, adjusted, color = symbol))      geom_point()      geom_line()      facet_wrap(~ symbol, ncol = 2, scales = "free_y")      scale_color_tq()      theme_tq()      labs(title = "After transformation: Easier to Understand")

澳门新浦京娱乐场网站 25

带有重复索引的时间序列

在某些应用场景中,可能会存在多个观测数据落在同一个时间点上的情况。下面就是一个例子:

In [63]: dates = pd.DatetimeIndex(['1/1/2000', '1/2/2000', '1/2/2000',
   ....:                           '1/2/2000', '1/3/2000'])
In [64]: dup_ts = pd.Series(np.arange(5), index=dates)

In [65]: dup_ts
Out[65]: 
2000-01-01    0
2000-01-02    1
2000-01-02    2
2000-01-02    3
2000-01-03    4
dtype: int64

通过检查索引的is_unique属性,我们就可以知道它是不是唯一的:

In [66]: dup_ts.index.is_unique
Out[66]: False

对这个时间序列进行索引,要么产生标量值,要么产生切片,具体要看所选的时间点是否重复:

In [67]: dup_ts['1/3/2000']  # not duplicated
Out[67]: 4

In [68]: dup_ts['1/2/2000']  # duplicated
Out[68]: 
2000-01-02    1
2000-01-02    2
2000-01-02    3
dtype: int64

假设你想要对具有非唯一时间戳的数据进行聚合。一个办法是使用groupby,并传入level=0:

In [69]: grouped = dup_ts.groupby(level=0)

In [70]: grouped.mean()
Out[70]: 
2000-01-01    0
2000-01-02    2
2000-01-03    4
dtype: int64

In [71]: grouped.count()
Out[71]: 
2000-01-01    1
2000-01-02    3
2000-01-03    1
dtype: int64

STEP 4:预测新数据

predict() 应用于回归模型。注意,和之前使用 lm() 函数时一样,去掉 indexdiff 列。

# Make predictionspred <- predict(    fit_lm,    newdata = select(        new_data_tbl, -c(index, diff)))predictions_tbl <- tibble(    date  = future_idx,    value = pred)predictions_tbl

## # A tibble: 12 x 2##          date     value##        <date>     <dbl>##  1 2017-01-01  9509.122##  2 2017-02-01 10909.189##  3 2017-03-01 12281.835##  4 2017-04-01 11378.678##  5 2017-05-01 13080.710##  6 2017-06-01 13583.471##  7 2017-07-01 11529.085##  8 2017-08-01 12984.939##  9 2017-09-01 11859.998## 10 2017-10-01 12331.419## 11 2017-11-01 13047.033## 12 2017-12-01 13940.003

sw_augument

sw_augument() 函数返回的 tibble 表中包含 .actual.fitted.resid 列,有助于在训练集上评估模型表现。注意,设置 timetk_idx = TRUE 返回初始的日期索引。

# sw_augment - get model residualssw_augment(fit_arima, timetk_idx = TRUE)

## # A tibble: 84 x 4##         index .actual   .fitted    .resid##        <date>   <dbl>     <dbl>     <dbl>##  1 2010-01-01    6558  6551.474  6.525878##  2 2010-02-01    7481  7473.583  7.416765##  3 2010-03-01    9475  9465.621  9.378648##  4 2010-04-01    9424  9414.704  9.295526##  5 2010-05-01    9351  9341.810  9.190414##  6 2010-06-01   10552 10541.641 10.359293##  7 2010-07-01    9077  9068.148  8.852178##  8 2010-08-01    9273  9263.984  9.016063##  9 2010-09-01    9420  9410.869  9.130943## 10 2010-10-01    9413  9403.908  9.091831## # ... with 74 more rows

我们可以可视化训练数据上的残差,看一下数据中有没有遗漏的模式没有被发现。

# Plotting residualssw_augment(fit_arima, timetk_idx = TRUE) %>%    ggplot(aes(x = index, y = .resid))      geom_point()       geom_hline(yintercept = 0, color = "red")       labs(title = "Residual diagnostic")      scale_x_date(date_breaks = "1 year", date_labels = "%Y")      theme_tq()

澳门新浦京娱乐场网站 26

tq_mutate

tq_mutate() 函数基于 xts 包为数据添加新的列。正因为这样,当返回数据不止一列时,tq_mutate() 显得特别有用(dplyr::mutate() 就没有这样的功能)。

11.3 日期的范围、频率以及移动

pandas中的原生时间序列一般被认为是不规则的,也就是说,它们没有固定的频率。对于大部分应用程序而言,这是无所谓的。但是,它常常需要以某种相对固定的频率进行分析,比如每日、每月、每15分钟等(这样自然会在时间序列中引入缺失值)。幸运的是,pandas有一整套标准时间序列频率以及用于重采样、频率推断、生成固定频率日期范围的工具。例如,我们可以将之前那个时间序列转换为一个具有固定频率(每日)的时间序列,只需调用resample即可:

In [72]: ts
Out[72]: 
2011-01-02   -0.204708
2011-01-05    0.478943
2011-01-07   -0.519439
2011-01-08   -0.555730
2011-01-10    1.965781
2011-01-12    1.393406
dtype: float64

In [73]: resampler = ts.resample('D')

字符串“D”是每天的意思。

频率的转换(或重采样)是一个比较大的主题,稍后将专门用一节来进行讨论(11.6小节)。这里,我将告诉你如何使用基本的频率和它的倍数。

STEP 5:比较实际值与预测值

我们可以使用 tq_get() 来检索实际数据。注意,我们没有用于比较的完整数据,但我们至少可以比较前几个月的实际值。

actuals_tbl <- tq_get(    "S4248SM144NCEN",    get = "economic.data",    from = "2017-01-01",    to = "2017-12-31")

可视化我们的预测。

# Plot Beer Sales Forecastbeer_sales_tbl %>%    ggplot(aes(x = date, y = price))      # Training data    geom_line(color = palette_light      geom_point(color = palette_light      # Predictions    geom_line(        aes(y = value),        color = palette_light()[[2]],        data = predictions_tbl)      geom_point(        aes(y = value),        color = palette_light()[[2]],        data = predictions_tbl)      # Actuals    geom_line(        color = palette_light()[[1]],        data = actuals_tbl)      geom_point(        color = palette_light()[[1]],        data = actuals_tbl)      # Aesthetics    theme_tq()      labs(        title = "Beer Sales Forecast: Time Series Machine Learning",        subtitle = "Using basic multivariate linear regression can yield accurate results")

澳门新浦京娱乐场网站 27

我们可以检查测试集上的错误(实际值 vs 预测值)。

# Investigate test errorerror_tbl <- left_join(    actuals_tbl, predictions_tbl) %>%    rename(        actual = price, pred = value) %>%    mutate(        error     = actual - pred,        error_pct = error / actual)         error_tbl

## # A tibble: 8 x 5##         date actual      pred     error    error_pct##       <date>  <int>     <dbl>     <dbl>        <dbl>## 1 2017-01-01   8664  9509.122 -845.1223 -0.097544127## 2 2017-02-01  10017 10909.189 -892.1891 -0.089067495## 3 2017-03-01  11960 12281.835 -321.8352 -0.026909301## 4 2017-04-01  11019 11378.678 -359.6777 -0.032641592## 5 2017-05-01  12971 13080.710 -109.7099 -0.008458092## 6 2017-06-01  14113 13583.471  529.5292  0.037520667## 7 2017-07-01  10928 11529.085 -601.0854 -0.055004156## 8 2017-08-01  12788 12984.939 -196.9386 -0.015400265

接着,我们可以做一些误差度量。预测值和实际值的 MAPE接近 4.5%,对于简单的多元线性回归模型来说是相当好的结果。更复杂的算法可以产生更精确的结果。

# Calculating test error metricstest_residuals <- error_tbl$errortest_error_pct <- error_tbl$error_pct * 100 # Percentage errorme   <- mean(test_residuals, na.rm=TRUE)rmse <- mean(test_residuals^2, na.rm=TRUE)^0.5mae  <- mean(abs(test_residuals), na.rm=TRUE)mape <- mean(abs(test_error_pct), na.rm=TRUE)mpe  <- mean(test_error_pct, na.rm=TRUE)tibble(me, rmse, mae, mape, mpe) %>% glimpse()

## Observations: 1## Variables: 5## $ me   <dbl> -349.6286## $ rmse <dbl> 551.7818## $ mae  <dbl> 482.0109## $ mape <dbl> 4.531821## $ mpe  <dbl> -3.593805

时间序列机器学习可能产生异常预测。对这个议题感兴趣的人可以阅读我们的短文 time series forecasting using timetk。

STEP 3:预测

使用 forecast() 函数做预测。

# Forecast next 12 monthsfcast_arima <- forecast(fit_arima, h = 12)

一个问题是,预测结果并不“tidy”。我们需要数据框形式的预测结果,以便应用 tidyverse 的功能,然而预测结果是 forecast 类型的,一种基于 ts 的对象。

class(fcast_arima)

## [1] "forecast"

tq_mutate 与滞后数据

一个关于 lag.xts 的例子。通常我们需要不只一列滞后数据,这正是 tq_mutate() 擅长的。下面,为原数据添加五列滞后数据。

# Lags - Get first 5 lags# Pro Tip: Make the new column names first, then add to the `col_rename` argcolumn_names <- paste0("lag", 1:5)# First five lags are output for each group of symbolsFANG_data_d %>%    select(symbol, date, adjusted) %>%    group_by %>%    tq_mutate(        select     = adjusted,        mutate_fun = lag.xts,        k          = 1:5,        col_rename = column_names)

## # A tibble: 3,024 x 8## # Groups:   symbol [4]##    symbol       date adjusted  lag1  lag2  lag3  lag4  lag5##     <chr>     <date>    <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>##  1     FB 2014-01-02    54.71    NA    NA    NA    NA    NA##  2     FB 2014-01-03    54.56 54.71    NA    NA    NA    NA##  3     FB 2014-01-06    57.20 54.56 54.71    NA    NA    NA##  4     FB 2014-01-07    57.92 57.20 54.56 54.71    NA    NA##  5     FB 2014-01-08    58.23 57.92 57.20 54.56 54.71    NA##  6     FB 2014-01-09    57.22 58.23 57.92 57.20 54.56 54.71##  7     FB 2014-01-10    57.94 57.22 58.23 57.92 57.20 54.56##  8     FB 2014-01-13    55.91 57.94 57.22 58.23 57.92 57.20##  9     FB 2014-01-14    57.74 55.91 57.94 57.22 58.23 57.92## 10     FB 2014-01-15    57.60 57.74 55.91 57.94 57.22 58.23## # ... with 3,014 more rows

生成日期范围

虽然我之前用的时候没有明说,但你可能已经猜到pandas.date_range可用于根据指定的频率生成指定长度的DatetimeIndex:

In [74]: index = pd.date_range('2012-04-01', '2012-06-01')

In [75]: index
Out[75]: 
DatetimeIndex(['2012-04-01', '2012-04-02', '2012-04-03', '2012-04-04',
               '2012-04-05', '2012-04-06', '2012-04-07', '2012-04-08',
               '2012-04-09', '2012-04-10', '2012-04-11', '2012-04-12',
               '2012-04-13', '2012-04-14', '2012-04-15', '2012-04-16',
               '2012-04-17', '2012-04-18', '2012-04-19', '2012-04-20',
               '2012-04-21', '2012-04-22', '2012-04-23', '2012-04-24',
               '2012-04-25', '2012-04-26', '2012-04-27', '2012-04-28',
               '2012-04-29', '2012-04-30', '2012-05-01', '2012-05-02',
               '2012-05-03', '2012-05-04', '2012-05-05', '2012-05-06',
               '2012-05-07', '2012-05-08', '2012-05-09', '2012-05-10',
               '2012-05-11', '2012-05-12', '2012-05-13', '2012-05-14',
               '2012-05-15', '2012-05-16', '2012-05-17', '2012-05-18',
               '2012-05-19', '2012-05-20', '2012-05-21', '2012-05-22',
               '2012-05-23', '2012-05-24', '2012-05-25', '2012-05-26',
               '2012-05-27', '2012-05-28', '2012-05-29', '2012-05-30',
               '2012-05-31', '2012-06-01'],
              dtype='datetime64[ns]', freq='D')

默认情况下,date_range会产生按天计算的时间点。如果只传入起始或结束日期,那就还得传入一个表示一段时间的数字:

In [76]: pd.date_range(start='2012-04-01', periods=20)
Out[76]: 
DatetimeIndex(['2012-04-01', '2012-04-02', '2012-04-03', '2012-04-04',
               '2012-04-05', '2012-04-06', '2012-04-07', '2012-04-08',
               '2012-04-09', '2012-04-10', '2012-04-11', '2012-04-12',
               '2012-04-13', '2012-04-14', '2012-04-15', '2012-04-16',
               '2012-04-17', '2012-04-18', '2012-04-19', '2012-04-20'],
              dtype='datetime64[ns]', freq='D')

In [77]: pd.date_range(end='2012-06-01', periods=20)
Out[77]: 
DatetimeIndex(['2012-05-13', '2012-05-14', '2012-05-15', '2012-05-16',
               '2012-05-17', '2012-05-18', '2012-05-19', '2012-05-20',
               '2012-05-21', '2012-05-22', '2012-05-23', '2012-05-24',
               '2012-05-25', '2012-05-26', '2012-05-27','2012-05-28',
               '2012-05-29', '2012-05-30', '2012-05-31', '2012-06-01'],
              dtype='datetime64[ns]', freq='D')

起始和结束日期定义了日期索引的严格边界。例如,如果你想要生成一个由每月最后一个工作日组成的日期索引,可以传入"BM"频率(表示business end of month,表11-4是频率列表),这样就只会包含时间间隔内(或刚好在边界上的)符合频率要求的日期:

In [78]: pd.date_range('2000-01-01', '2000-12-01', freq='BM')
Out[78]: 
DatetimeIndex(['2000-01-31', '2000-02-29', '2000-03-31', '2000-04-28',
               '2000-05-31', '2000-06-30', '2000-07-31', '2000-08-31',
               '2000-09-29', '2000-10-31', '2000-11-30'],
              dtype='datetime64[ns]', freq='BM')

表11-4 基本的时间序列频率(不完整)

澳门新浦京娱乐场网站 28

澳门新浦京娱乐场网站 29

澳门新浦京娱乐场网站 30

date_range默认会保留起始和结束时间戳的时间信息(如果有的话):

In [79]: pd.date_range('2012-05-02 12:56:31', periods=5)
Out[79]: 
DatetimeIndex(['2012-05-02 12:56:31', '2012-05-03 12:56:31',
               '2012-05-04 12:56:31', '2012-05-05 12:56:31',
               '2012-05-06 12:56:31'],
              dtype='datetime64[ns]', freq='D')

有时,虽然起始和结束日期带有时间信息,但你希望产生一组被规范化(normalize)到午夜的时间戳。normalize选项即可实现该功能:

In [80]: pd.date_range('2012-05-02 12:56:31', periods=5, normalize=True)
Out[80]: 
DatetimeIndex(['2012-05-02', '2012-05-03', '2012-05-04', '2012-05-05',
               '2012-05-06'],
              dtype='datetime64[ns]', freq='D')

PART 2:转换

  • 问题:R 中不同类型的时间序列数据难以方便一致地实现相互转换。
  • 解决方案tk_tbltk_xtstk_zootk_ts

STEP 4:用 sweep 简化预测

我们使用 sw_sweep() 简化预测结果,一个额外的好处是,如果 forecast 对象有 timetk 索引,我们可以用它返回一个日期时间索引,不同于 ts 对象的规则索引。

首先要确认 forecast 对象有 timetk 索引,这需要在使用 sw_sweep() 时设置 timetk_idx 参数。

# Check if object has timetk index has_timetk_idx(fcast_arima)

## [1] TRUE

现在,使用 sw_sweep() 来简化预测结果,它会在内部根据 time_tk 构造一条未来时间序列索引(这一步总是会被执行,因为我们在第 1 步中用 tk_ts() 构造了 ts 对象)注意:这意味着我们最终可以在 forecast 包中使用日期(不同于 ts 对象中的规则索引)!

# sw_sweep - tidies forecast outputfcast_tbl <- sw_sweep(fcast_arima, timetk_idx = TRUE)fcast_tbl

## # A tibble: 96 x 7##         index    key price lo.80 lo.95 hi.80 hi.95##        <date>  <chr> <dbl> <dbl> <dbl> <dbl> <dbl>##  1 2010-01-01 actual  6558    NA    NA    NA    NA##  2 2010-02-01 actual  7481    NA    NA    NA    NA##  3 2010-03-01 actual  9475    NA    NA    NA    NA##  4 2010-04-01 actual  9424    NA    NA    NA    NA##  5 2010-05-01 actual  9351    NA    NA    NA    NA##  6 2010-06-01 actual 10552    NA    NA    NA    NA##  7 2010-07-01 actual  9077    NA    NA    NA    NA##  8 2010-08-01 actual  9273    NA    NA    NA    NA##  9 2010-09-01 actual  9420    NA    NA    NA    NA## 10 2010-10-01 actual  9413    NA    NA    NA    NA## # ... with 86 more rows

tq_mutate 与滚动函数

另一个例子,应用 xts 中的滚动函数 roll.apply()。让我们借助函数 quantile() 得到滚动分位数。下面是每个函数的参数:

tq_mutate 的参数:

  • select = adjusted,只选择复权修正过的数据列。该参数也可以不填,或选择其他不同的列。
  • mutate_fun = rollapply,这是一个 xts 函数,将会以 “tidy” 的方式调用。

rollapply 的参数:

  • width = 5,告诉 rollapply 计算窗口的周期是多少。
  • by.column = FALSErollapply() 函数默认对每一列分别操作,然而我们要把所有列放在一起操作。
  • FUN = quantilequantile() 正是要被滚动调用的函数。

quantile 的参数:

  • probs = c(0, 0.025, ...),计算这些概率的分位数。
  • na.rm = TRUEquantile 会去掉遇到的 NA 值。
# Rolling quantileFANG_data_d %>%    select(symbol, date, adjusted) %>%    group_by %>%    tq_mutate(        select     = adjusted,        mutate_fun = rollapply,        width      = 5,        by.column  = FALSE,        FUN        = quantile,        probs      = c(0, 0.025, 0.25, 0.5, 0.75, 0.975, 1),        na.rm      = TRUE)

## # A tibble: 3,024 x 10## # Groups:   symbol [4]##    symbol       date adjusted   X0.  X2.5.  X25.  X50.  X75. X97.5.##     <chr>     <date>    <dbl> <dbl>  <dbl> <dbl> <dbl> <dbl>  <dbl>##  1     FB 2014-01-02    54.71    NA     NA    NA    NA    NA     NA##  2     FB 2014-01-03    54.56    NA     NA    NA    NA    NA     NA##  3     FB 2014-01-06    57.20    NA     NA    NA    NA    NA     NA##  4     FB 2014-01-07    57.92    NA     NA    NA    NA    NA     NA##  5     FB 2014-01-08    58.23 54.56 54.575 54.71 57.20 57.92 58.199##  6     FB 2014-01-09    57.22 54.56 54.824 57.20 57.22 57.92 58.199##  7     FB 2014-01-10    57.94 57.20 57.202 57.22 57.92 57.94 58.201##  8     FB 2014-01-13    55.91 55.91 56.041 57.22 57.92 57.94 58.201##  9     FB 2014-01-14    57.74 55.91 56.041 57.22 57.74 57.94 58.201## 10     FB 2014-01-15    57.60 55.91 56.041 57.22 57.60 57.74 57.920## # ... with 3,014 more rows, and 1 more variables: X100. <dbl>

频率和日期偏移量

pandas中的频率是由一个基础频率(base frequency)和一个乘数组成的。基础频率通常以一个字符串别名表示,比如"M"表示每月,"H"表示每小时。对于每个基础频率,都有一个被称为日期偏移量(date offset)的对象与之对应。例如,按小时计算的频率可以用Hour类表示:

In [81]: from pandas.tseries.offsets import Hour, Minute

In [82]: hour = Hour()

In [83]: hour
Out[83]: <Hour>

传入一个整数即可定义偏移量的倍数:

In [84]: four_hours = Hour(4)

In [85]: four_hours
Out[85]: <4 * Hours>

一般来说,无需明确创建这样的对象,只需使用诸如"H"或"4H"这样的字符串别名即可。在基础频率前面放上一个整数即可创建倍数:

In [86]: pd.date_range('2000-01-01', '2000-01-03 23:59', freq='4h')
Out[86]: 
DatetimeIndex(['2000-01-01 00:00:00', '2000-01-01 04:00:00',
               '2000-01-01 08:00:00', '2000-01-01 12:00:00',
               '2000-01-01 16:00:00', '2000-01-01 20:00:00',
               '2000-01-02 00:00:00', '2000-01-02 04:00:00',
               '2000-01-02 08:00:00', '2000-01-02 12:00:00',
               '2000-01-02 16:00:00', '2000-01-02 20:00:00',
               '2000-01-03 00:00:00', '2000-01-03 04:00:00',
               '2000-01-03 08:00:00', '2000-01-03 12:00:00',
               '2000-01-03 16:00:00', '2000-01-03 20:00:00'],
              dtype='datetime64[ns]', freq='4H')

大部分偏移量对象都可通过加法进行连接:

In [87]: Hour(2)   Minute(30)
Out[87]: <150 * Minutes>

同理,你也可以传入频率字符串(如"2h30min"),这种字符串可以被高效地解析为等效的表达式:

In [88]: pd.date_range('2000-01-01', periods=10, freq='1h30min')
Out[88]: 
DatetimeIndex(['2000-01-01 00:00:00', '2000-01-01 01:30:00',
               '2000-01-01 03:00:00', '2000-01-01 04:30:00',
               '2000-01-01 06:00:00', '2000-01-01 07:30:00',
               '2000-01-01 09:00:00', '2000-01-01 10:30:00',
               '2000-01-01 12:00:00', '2000-01-01 13:30:00'],
              dtype='datetime64[ns]', freq='90T')

有些频率所描述的时间点并不是均匀分隔的。例如,"M"(日历月末)和"BM"(每月最后一个工作日)就取决于每月的天数,对于后者,还要考虑月末是不是周末。由于没有更好的术语,我将这些称为锚点偏移量(anchored offset)。

表11-4列出了pandas中的频率代码和日期偏移量类。

笔记:用户可以根据实际需求自定义一些频率类以便提供pandas所没有的日期逻辑,但具体的细节超出了本书的范围。

表11-4 时间序列的基础频率

澳门新浦京娱乐场网站 31

澳门新浦京娱乐场网站 32

澳门新浦京娱乐场网站 33

tk_xts

我们开始的时候用 tbl 对象,一个劣势是有时候必须转换成 xts 对象,因为要使用其他包(例如 xtszooquantmod 等等)里基于 xts 对象的函数。

我们可以使用 tk_xts() 函数轻松地将数据转换成 xts 对象。注意,tk_xts() 函数会自动检测包含时间的列,并把该列当做 xts 对象的索引。

# Coerce to xtsbeer_sales_xts <- tk_xts(beer_sales_tbl) # Show the first six rows of the xts objectbeer_sales_xts %>%    head()

##            price## 2010-01-01  6558## 2010-02-01  7481## 2010-03-01  9475## 2010-04-01  9424## 2010-05-01  9351## 2010-06-01 10552

我们也可以从 xts 转成 tbl 。我们设定 rename_index = "date" 让索引的名字和开始的时候保持一致。这种操作在以前不太容易。

tk_tbl(beer_sales_xts, rename_index = "date")

## # A tibble: 84 x 2##          date price##        <date> <int>##  1 2010-01-01  6558##  2 2010-02-01  7481##  3 2010-03-01  9475##  4 2010-04-01  9424##  5 2010-05-01  9351##  6 2010-06-01 10552##  7 2010-07-01  9077##  8 2010-08-01  9273##  9 2010-09-01  9420## 10 2010-10-01  9413## # ... with 74 more rows

STEP 5:比较真实值和预测值

我们可以使用 tq_get() 来检索实际数据。注意,我们没有用于比较的完整数据,但我们至少可以比较前几个月的实际值。

actuals_tbl <- tq_get(    "S4248SM144NCEN",    get = "economic.data",    from = "2017-01-01",    to = "2017-12-31")

注意,预测结果放在 tibble 中,可以方便的实现可视化。

# Visualize the forecast with ggplotfcast_tbl %>%    ggplot(        aes(x = index, y = price, color = key))      # 95% CI    geom_ribbon(        aes(ymin = lo.95, ymax = hi.95),         fill = "#D5DBFF", color = NA, size = 0)      # 80% CI    geom_ribbon(        aes(ymin = lo.80, ymax = hi.80, fill = key),         fill = "#596DD5", color = NA,         size = 0, alpha = 0.8)      # Prediction    geom_line()      geom_point()      # Actuals    geom_line(        aes(x = date, y = price), color = palette_light()[[1]],        data = actuals_tbl)      geom_point(        aes(x = date, y = price), color = palette_light()[[1]],         data = actuals_tbl)      # Aesthetics    labs(        title = "Beer Sales Forecast: ARIMA", x = "", y = "Thousands of Tons",        subtitle = "sw_sweep tidies the auto.arima() forecast output")      scale_x_date(        date_breaks = "1 year",        date_labels = "%Y")      scale_color_tq()      scale_fill_tq()      theme_tq()

澳门新浦京娱乐场网站 34

我们可以研究测试集上的误差(真实值 vs 预测值)。

# Investigate test error error_tbl <- left_join(    actuals_tbl,    fcast_tbl,    by = c("date" = "index")) %>%    rename(        actual = price.x, pred = price.y) %>%    select(date, actual, pred) %>%    mutate(        error     = actual - pred,        error_pct = error / actual)         error_tbl

## # A tibble: 8 x 5##         date actual      pred      error    error_pct##       <date>  <int>     <dbl>      <dbl>        <dbl>## 1 2017-01-01   8664  8601.815   62.18469  0.007177365## 2 2017-02-01  10017 10855.429 -838.42908 -0.083700617## 3 2017-03-01  11960 11502.214  457.78622  0.038276439## 4 2017-04-01  11019 11582.600 -563.59962 -0.051147982## 5 2017-05-01  12971 12566.765  404.23491  0.031164514## 6 2017-06-01  14113 13263.918  849.08191  0.060163106## 7 2017-07-01  10928 11507.277 -579.27693 -0.053008504## 8 2017-08-01  12788 12527.278  260.72219  0.020388035

并且,我们可以做简单的误差度量。MAPE 接近 4.3%,比简单的线性回归模型略好一点,但是 RMSE 变差了。

# Calculate test error metricstest_residuals <- error_tbl$errortest_error_pct <- error_tbl$error_pct * 100 # Percentage errorme   <- mean(test_residuals, na.rm=TRUE)rmse <- mean(test_residuals^2, na.rm=TRUE)^0.5mae  <- mean(abs(test_residuals), na.rm=TRUE)mape <- mean(abs(test_error_pct), na.rm=TRUE)mpe  <- mean(test_error_pct, na.rm=TRUE)tibble(me, rmse, mae, mape, mpe) %>%    glimpse()

## Observations: 1## Variables: 5## $ me   <dbl> 6.588034## $ rmse <dbl> 561.4631## $ mae  <dbl> 501.9144## $ mape <dbl> 4.312832## $ mpe  <dbl> -0.3835956

可用函数

已经介绍了如何将 xts 函数和 tq_transmutetq_mutate 联合使用。还有许多 xts 函数可以以 “tidy” 的方式使用!用 tq_transmute_fun_options() 查看其他可用函数。

# Available functions# mutate_fun =tq_transmute_fun_options()

## $zoo##  [1] "rollapply"          "rollapplyr"         "rollmax"           ##  [4] "rollmax.default"    "rollmaxr"           "rollmean"          ##  [7] "rollmean.default"   "rollmeanr"          "rollmedian"        ## [10] "rollmedian.default" "rollmedianr"        "rollsum"           ## [13] "rollsum.default"    "rollsumr"          ## ## $xts##  [1] "apply.daily"     "apply.monthly"   "apply.quarterly"##  [4] "apply.weekly"    "apply.yearly"    "diff.xts"       ##  [7] "lag.xts"         "period.apply"    "period.max"     ## [10] "period.min"      "period.prod"     "period.sum"     ## [13] "periodicity"     "to.daily"        "to.hourly"      ## [16] "to.minutes"      "to.minutes10"    "to.minutes15"   ## [19] "to.minutes3"     "to.minutes30"    "to.minutes5"    ## [22] "to.monthly"      "to.period"       "to.quarterly"   ## [25] "to.weekly"       "to.yearly"       "to_period"      ## ## $quantmod##  [1] "allReturns"      "annualReturn"    "ClCl"           ##  [4] "dailyReturn"     "Delt"            "HiCl"           ##  [7] "Lag"             "LoCl"            "LoHi"           ## [10] "monthlyReturn"   "Next"            "OpCl"           ## [13] "OpHi"            "OpLo"            "OpOp"           ## [16] "periodReturn"    "quarterlyReturn" "seriesAccel"    ## [19] "seriesDecel"     "seriesDecr"      "seriesHi"       ## [22] "seriesIncr"      "seriesLo"        "weeklyReturn"   ## [25] "yearlyReturn"   ## ## $TTR##  [1] "adjRatios"          "ADX"                "ALMA"              ##  [4] "aroon"              "ATR"                "BBands"            ##  [7] "CCI"                "chaikinAD"          "chaikinVolatility" ## [10] "CLV"                "CMF"                "CMO"               ## [13] "DEMA"               "DonchianChannel"    "DPO"               ## [16] "DVI"                "EMA"                "EMV"               ## [19] "EVWMA"              "GMMA"               "growth"            ## [22] "HMA"                "KST"                "lags"              ## [25] "MACD"               "MFI"                "momentum"          ## [28] "OBV"                "PBands"             "ROC"               ## [31] "rollSFM"            "RSI"                "runCor"            ## [34] "runCov"             "runMAD"             "runMax"            ## [37] "runMean"            "runMedian"          "runMin"            ## [40] "runPercentRank"     "runSD"              "runSum"            ## [43] "runVar"             "SAR"                "SMA"               ## [46] "SMI"                "SNR"                "stoch"             ## [49] "TDI"                "TRIX"               "ultimateOscillator"## [52] "VHF"                "VMA"                "volatility"        ## [55] "VWAP"               "VWMA"               "wilderSum"         ## [58] "williamsAD"         "WMA"                "WPR"               ## [61] "ZigZag"             "ZLEMA"             ## ## $PerformanceAnalytics## [1] "Return.annualized"        "Return.annualized.excess"## [3] "Return.clean"             "Return.cumulative"       ## [5] "Return.excess"            "Return.Geltner"          ## [7] "zerofill"

WOM日期

WOM(Week Of Month)是一种非常实用的频率类,它以WOM开头。它使你能获得诸如“每月第3个星期五”之类的日期:

In [89]: rng = pd.date_range('2012-01-01', '2012-09-01', freq='WOM-3FRI')

In [90]: list(rng)
Out[90]: 
[Timestamp('2012-01-20 00:00:00', freq='WOM-3FRI'),
 Timestamp('2012-02-17 00:00:00', freq='WOM-3FRI'),
 Timestamp('2012-03-16 00:00:00', freq='WOM-3FRI'),
 Timestamp('2012-04-20 00:00:00', freq='WOM-3FRI'),
 Timestamp('2012-05-18 00:00:00', freq='WOM-3FRI'),
 Timestamp('2012-06-15 00:00:00', freq='WOM-3FRI'),
 Timestamp('2012-07-20 00:00:00', freq='WOM-3FRI'),
 Timestamp('2012-08-17 00:00:00', freq='WOM-3FRI')]

tk_ts

有许多包用了另一种类型的时间序列数据——ts,其中最常见的可能就是 forecast 包。使用 tk_ts() 函数的优点有两个:

  1. 与其他 tk_ 函数兼容,可以方便直接的实现数据的转换和逆转换。
  2. 更重要的是:当使用 tk_ts 函数时,ts 对象将初始的不规则时间索引转换成一个索引属性。这使得保留日期和时间信息成为可能。

一个例子,使用 tk_tbl(timetk_index = TRUE)函数转换成 ts 对象。因为 ts 对象只能处理规则时间索引,我们需要添加参数 start = 2010freq = 12

# Coerce to tsbeer_sales_ts <- tk_ts(    beer_sales_tbl,    start = 2010,    freq = 12)# Show the calendar-printoutbeer_sales_ts

##        Jan   Feb   Mar   Apr   May   Jun   Jul   Aug   Sep   Oct## 2010  6558  7481  9475  9424  9351 10552  9077  9273  9420  9413## 2011  6901  8014  9833  9281  9967 11344  9106 10468 10085  9612## 2012  7486  8641  9709  9423 11342 11274  9845 11163  9532 10754## 2013  8395  8888 10109 10493 12217 11385 11186 11462 10494 11541## 2014  8559  9061 10058 10979 11794 11906 10966 10981 10827 11815## 2015  8398  9061 10720 11105 11505 12903 11866 11223 12023 11986## 2016  8540 10158 11879 11155 11916 13291 10540 12212 11786 11424##        Nov   Dec## 2010  9866 11455## 2011 10328 11483## 2012 10953 11922## 2013 11139 12709## 2014 10466 13303## 2015 11510 14190## 2016 12482 13832

有两种方法转换回 tbl

  1. 直接使用 tk_tbl(),我们将得到 YEARMON 类型的“规则”的时间索引(来自 zoo 包)。
  2. 如果原始对象用 tk_ts() 创建,并且有属性 timetk_index,我们可以通过命令 tk_tbl(timetk_index = TRUE) 转换回去,并得到 Date 格式 “非规则”时间索引。

方法 1:注意,日期列是 YEARMON 类型的。

tk_tbl(beer_sales_ts, rename_index = "date")

## # A tibble: 84 x 2##             date price##    <S3: yearmon> <int>##  1      Jan 2010  6558##  2      Feb 2010  7481##  3      Mar 2010  9475##  4      Apr 2010  9424##  5      May 2010  9351##  6      Jun 2010 10552##  7      Jul 2010  9077##  8      Aug 2010  9273##  9      Sep 2010  9420## 10      Oct 2010  9413## # ... with 74 more rows

方法 2:设置参数 timetk_idx = TRUE找回初始的日期或时间信息

首先,用 has_timetk_idx() 检查 ts 对象是否存在 timetk 索引。

# Check for timetk index. has_timetk_idx(beer_sales_ts)

## [1] TRUE

如果返回值是 TRUE,在调用 tk_tbl() 时设定 timetk_idx = TRUE。现在可以看到日期列是 date 类型的,这在以往不容易做到。

# If timetk_idx is present, can get original dates back tk_tbl(beer_sales_ts, timetk_idx = TRUE, rename_index = "date")

## # A tibble: 84 x 2##          date price##        <date> <int>##  1 2010-01-01  6558##  2 2010-02-01  7481##  3 2010-03-01  9475##  4 2010-04-01  9424##  5 2010-05-01  9351##  6 2010-06-01 10552##  7 2010-07-01  9077##  8 2010-08-01  9273##  9 2010-09-01  9420## 10 2010-10-01  9413## # ... with 74 more rows

移动(超前和滞后)数据

移动(shifting)指的是沿着时间轴将数据前移或后移。Series和DataFrame都有一个shift方法用于执行单纯的前移或后移操作,保持索引不变:

In [91]: ts = pd.Series(np.random.randn(4),
   ....:                index=pd.date_range('1/1/2000', periods=4, freq='M'))

In [92]: ts
Out[92]: 
2000-01-31   -0.066748
2000-02-29    0.838639
2000-03-31   -0.117388
2000-04-30   -0.517795
Freq: M, dtype: float64

In [93]: ts.shift(2)
Out[93]: 
2000-01-31         NaN
2000-02-29         NaN
2000-03-31   -0.066748
2000-04-30    0.838639
Freq: M, dtype: float64

In [94]: ts.shift(-2)
Out[94]: 
2000-01-31   -0.117388
2000-02-29   -0.517795
2000-03-31         NaN
2000-04-30         NaN
Freq: M, dtype: float64

当我们这样进行移动时,就会在时间序列的前面或后面产生缺失数据。

shift通常用于计算一个时间序列或多个时间序列(如DataFrame的列)中的百分比变化。可以这样表达:

ts / ts.shift(1) - 1

由于单纯的移位操作不会修改索引,所以部分数据会被丢弃。因此,如果频率已知,则可以将其传给shift以便实现对时间戳进行位移而不是对数据进行简单位移:

In [95]: ts.shift(2, freq='M')
Out[95]: 
2000-03-31   -0.066748
2000-04-30    0.838639
2000-05-31   -0.117388
2000-06-30   -0.517795
Freq: M, dtype: float64

这里还可以使用其他频率,于是你就能非常灵活地对数据进行超前和滞后处理了:

In [96]: ts.shift(3, freq='D')
Out[96]: 
2000-02-03   -0.066748
2000-03-03    0.838639
2000-04-03   -0.117388
2000-05-03   -0.517795
dtype: float64

In [97]: ts.shift(1, freq='90T')
Out[97]: 
2000-01-31 01:30:00   -0.066748
2000-02-29 01:30:00    0.838639
2000-03-31 01:30:00   -0.117388
2000-04-30 01:30:00   -0.517795
Freq: M, dtype: float64

通过偏移量对日期进行位移

pandas的日期偏移量还可以用在datetime或Timestamp对象上:

In [98]: from pandas.tseries.offsets import Day, MonthEnd

In [99]: now = datetime(2011, 11, 17)

In [100]: now   3 * Day()
Out[100]: Timestamp('2011-11-20 00:00:00')

如果加的是锚点偏移量(比如MonthEnd),第一次增量会将原日期向前滚动到符合频率规则的下一个日期:

In [101]: now   MonthEnd()
Out[101]: Timestamp('2011-11-30 00:00:00')

In [102]: now   MonthEnd(2)
Out[102]: Timestamp('2011-12-31 00:00:00')

通过锚点偏移量的rollforward和rollback方法,可明确地将日期向前或向后“滚动”:

In [103]: offset = MonthEnd()

In [104]: offset.rollforward(now)
Out[104]: Timestamp('2011-11-30 00:00:00')

In [105]: offset.rollback(now)
Out[105]: Timestamp('2011-10-31 00:00:00')

日期偏移量还有一个巧妙的用法,即结合groupby使用这两个“滚动”方法:

In [106]: ts = pd.Series(np.random.randn(20),
   .....:                index=pd.date_range('1/15/2000', periods=20, freq='4d'))

In [107]: ts
Out[107]: 
2000-01-15   -0.116696
2000-01-19    2.389645
2000-01-23   -0.932454
2000-01-27   -0.229331
2000-01-31   -1.140330
2000-02-04    0.439920
2000-02-08   -0.823758
2000-02-12   -0.520930
2000-02-16    0.350282
2000-02-20    0.204395
2000-02-24    0.133445
2000-02-28    0.327905
2000-03-03    0.072153
2000-03-07    0.131678
2000-03-11   -1.297459
2000-03-15    0.997747
2000-03-19    0.870955
2000-03-23   -0.991253
2000-03-27    0.151699
2000-03-31    1.266151
Freq: 4D, dtype: float64

In [108]: ts.groupby(offset.rollforward).mean()
Out[108]: 
2000-01-31   -0.005833
2000-02-29    0.015894
2000-03-31    0.150209
dtype: float64

当然,更简单、更快速地实现该功能的办法是使用resample(11.6小节将对此进行详细介绍):

In [109]: ts.resample('M').mean()
Out[109]: 
2000-01-31   -0.005833
2000-02-29    0.015894
2000-03-31    0.150209
Freq: M, dtype: float64

11.4 时区处理

时间序列处理工作中最让人不爽的就是对时区的处理。许多人都选择以协调世界时(UTC,它是格林尼治标准时间(Greenwich Mean Time)的接替者,目前已经是国际标准了)来处理时间序列。时区是以UTC偏移量的形式表示的。例如,夏令时期间,纽约比UTC慢4小时,而在全年其他时间则比UTC慢5小时。

在Python中,时区信息来自第三方库pytz,它使Python可以使用Olson数据库(汇编了世界时区信息)。这对历史数据非常重要,这是因为由于各地政府的各种突发奇想,夏令时转变日期(甚至UTC偏移量)已经发生过多次改变了。就拿美国来说,DST转变时间自1900年以来就改变过多次!

有关pytz库的更多信息,请查阅其文档。就本书而言,由于pandas包装了pytz的功能,因此你可以不用记忆其API,只要记得时区的名称即可。时区名可以在shell中看到,也可以通过文档查看:

In [110]: import pytz

In [111]: pytz.common_timezones[-5:]
Out[111]: ['US/Eastern', 'US/Hawaii', 'US/Mountain', 'US/Pacific', 'UTC']

要从pytz中获取时区对象,使用pytz.timezone即可:

In [112]: tz = pytz.timezone('America/New_York')

In [113]: tz
Out[113]: <DstTzInfo 'America/New_York' LMT-1 day, 19:04:00 STD>

pandas中的方法既可以接受时区名也可以接受这些对象。

时区本地化和转换

默认情况下,pandas中的时间序列是单纯的(naive)时区。看看下面这个时间序列:

In [114]: rng = pd.date_range('3/9/2012 9:30', periods=6, freq='D')

In [115]: ts = pd.Series(np.random.randn(len(rng)), index=rng)

In [116]: ts
Out[116]: 
2012-03-09 09:30:00   -0.202469
2012-03-10 09:30:00    0.050718
2012-03-11 09:30:00    0.639869
2012-03-12 09:30:00    0.597594
2012-03-13 09:30:00   -0.797246
2012-03-14 09:30:00    0.472879
Freq: D, dtype: float64

其索引的tz字段为None:

In [117]: print(ts.index.tz)
None

可以用时区集生成日期范围:

In [118]: pd.date_range('3/9/2012 9:30', periods=10, freq='D', tz='UTC')
Out[118]: 
DatetimeIndex(['2012-03-09 09:30:00 00:00', '2012-03-10 09:30:00 00:00',
               '2012-03-11 09:30:00 00:00', '2012-03-12 09:30:00 00:00',
               '2012-03-13 09:30:00 00:00', '2012-03-14 09:30:00 00:00',
               '2012-03-15 09:30:00 00:00', '2012-03-16 09:30:00 00:00',
               '2012-03-17 09:30:00 00:00', '2012-03-18 09:30:00 00:00'],
              dtype='datetime64[ns, UTC]', freq='D')

从单纯到本地化的转换是通过tz_localize方法处理的:

In [119]: ts
Out[119]: 
2012-03-09 09:30:00   -0.202469
2012-03-10 09:30:00    0.050718
2012-03-11 09:30:00    0.639869
2012-03-12 09:30:00    0.597594
2012-03-13 09:30:00   -0.797246
2012-03-14 09:30:00    0.472879
Freq: D, dtype: float64

In [120]: ts_utc = ts.tz_localize('UTC')

In [121]: ts_utc
Out[121]: 
2012-03-09 09:30:00 00:00   -0.202469
2012-03-10 09:30:00 00:00    0.050718
2012-03-11 09:30:00 00:00    0.639869
2012-03-12 09:30:00 00:00    0.597594
2012-03-13 09:30:00 00:00   -0.797246
2012-03-14 09:30:00 00:00    0.472879
Freq: D, dtype: float64

In [122]: ts_utc.index
Out[122]: 
DatetimeIndex(['2012-03-09 09:30:00 00:00', '2012-03-10 09:30:00 00:00',
               '2012-03-11 09:30:00 00:00', '2012-03-12 09:30:00 00:00',
               '2012-03-13 09:30:00 00:00', '2012-03-14 09:30:00 00:00'],
              dtype='datetime64[ns, UTC]', freq='D')

一旦时间序列被本地化到某个特定时区,就可以用tz_convert将其转换到别的时区了:

In [123]: ts_utc.tz_convert('America/New_York')
Out[123]: 
2012-03-09 04:30:00-05:00   -0.202469
2012-03-10 04:30:00-05:00    0.050718
2012-03-11 05:30:00-04:00    0.639869
2012-03-12 05:30:00-04:00    0.597594
2012-03-13 05:30:00-04:00   -0.797246
2012-03-14 05:30:00-04:00    0.472879
Freq: D, dtype: float64

对于上面这种时间序列(它跨越了美国东部时区的夏令时转变期),我们可以将其本地化到EST,然后转换为UTC或柏林时间:

In [124]: ts_eastern = ts.tz_localize('America/New_York')

In [125]: ts_eastern.tz_convert('UTC')
Out[125]: 
2012-03-09 14:30:00 00:00   -0.202469
2012-03-10 14:30:00 00:00    0.050718
2012-03-11 13:30:00 00:00    0.639869
2012-03-12 13:30:00 00:00    0.597594
2012-03-13 13:30:00 00:00   -0.797246
2012-03-14 13:30:00 00:00    0.472879
Freq: D, dtype: float64

In [126]: ts_eastern.tz_convert('Europe/Berlin')
Out[126]: 
2012-03-09 15:30:00 01:00   -0.202469
2012-03-10 15:30:00 01:00    0.050718
2012-03-11 14:30:00 01:00    0.639869
2012-03-12 14:30:00 01:00    0.597594
2012-03-13 14:30:00 01:00   -0.797246
2012-03-14 14:30:00 01:00    0.472879
Freq: D, dtype: float64

tz_localize和tz_convert也是DatetimeIndex的实例方法:

In [127]: ts.index.tz_localize('Asia/Shanghai')
Out[127]: 
DatetimeIndex(['2012-03-09 09:30:00 08:00', '2012-03-10 09:30:00 08:00',
               '2012-03-11 09:30:00 08:00', '2012-03-12 09:30:00 08:00',
               '2012-03-13 09:30:00 08:00', '2012-03-14 09:30:00 08:00'],
              dtype='datetime64[ns, Asia/Shanghai]', freq='D')

注意:对单纯时间戳的本地化操作还会检查夏令时转变期附近容易混淆或不存在的时间。

操作时区意识型Timestamp对象

跟时间序列和日期范围差不多,独立的Timestamp对象也能被从单纯型(naive)本地化为时区意识型(time zone-aware),并从一个时区转换到另一个时区:

In [128]: stamp = pd.Timestamp('2011-03-12 04:00')

In [129]: stamp_utc = stamp.tz_localize('utc')

In [130]: stamp_utc.tz_convert('America/New_York')
Out[130]: Timestamp('2011-03-11 23:00:00-0500', tz='America/New_York')

在创建Timestamp时,还可以传入一个时区信息:

In [131]: stamp_moscow = pd.Timestamp('2011-03-12 04:00', tz='Europe/Moscow')

In [132]: stamp_moscow
Out[132]: Timestamp('2011-03-12 04:00:00 0300', tz='Europe/Moscow')

时区意识型Timestamp对象在内部保存了一个UTC时间戳值(自UNIX纪元(1970年1月1日)算起的纳秒数)。这个UTC值在时区转换过程中是不会发生变化的:

In [133]: stamp_utc.value
Out[133]: 1299902400000000000

In [134]: stamp_utc.tz_convert('America/New_York').value
Out[134]: 1299902400000000000

当使用pandas的DateOffset对象执行时间算术运算时,运算过程会自动关注是否存在夏令时转变期。这里,我们创建了在DST转变之前的时间戳。首先,来看夏令时转变前的30分钟:

In [135]: from pandas.tseries.offsets import Hour

In [136]: stamp = pd.Timestamp('2012-03-12 01:30', tz='US/Eastern')

In [137]: stamp
Out[137]: Timestamp('2012-03-12 01:30:00-0400', tz='US/Eastern')

In [138]: stamp   Hour()
Out[138]: Timestamp('2012-03-12 02:30:00-0400', tz='US/Eastern')

然后,夏令时转变前90分钟:

In [139]: stamp = pd.Timestamp('2012-11-04 00:30', tz='US/Eastern')

In [140]: stamp
Out[140]: Timestamp('2012-11-04 00:30:00-0400', tz='US/Eastern')

In [141]: stamp   2 * Hour()
Out[141]: Timestamp('2012-11-04 01:30:00-0500', tz='US/Eastern')

不同时区之间的运算

如果两个时间序列的时区不同,在将它们合并到一起时,最终结果就会是UTC。由于时间戳其实是以UTC存储的,所以这是一个很简单的运算,并不需要发生任何转换:

In [142]: rng = pd.date_range('3/7/2012 9:30', periods=10, freq='B')

In [143]: ts = pd.Series(np.random.randn(len(rng)), index=rng)

In [144]: ts
Out[144]: 
2012-03-07 09:30:00    0.522356
2012-03-08 09:30:00   -0.546348
2012-03-09 09:30:00   -0.733537
2012-03-12 09:30:00    1.302736
2012-03-13 09:30:00    0.022199
2012-03-14 09:30:00    0.364287
2012-03-15 09:30:00   -0.922839
2012-03-16 09:30:00    0.312656
2012-03-19 09:30:00   -1.128497
2012-03-20 09:30:00   -0.333488
Freq: B, dtype: float64

In [145]: ts1 = ts[:7].tz_localize('Europe/London')

In [146]: ts2 = ts1[2:].tz_convert('Europe/Moscow')

In [147]: result = ts1   ts2

In [148]: result.index
Out[148]: 
DatetimeIndex(['2012-03-07 09:30:00 00:00', '2012-03-08 09:30:00 00:00',
               '2012-03-09 09:30:00 00:00', '2012-03-12 09:30:00 00:00',
               '2012-03-13 09:30:00 00:00', '2012-03-14 09:30:00 00:00',
               '2012-03-15 09:30:00 00:00'],
              dtype='datetime64[ns, UTC]', freq='B')

11.5 时期及其算术运算

时期(period)表示的是时间区间,比如数日、数月、数季、数年等。Period类所表示的就是这种数据类型,其构造函数需要用到一个字符串或整数,以及表11-4中的频率:

In [149]: p = pd.Period(2007, freq='A-DEC')

In [150]: p
Out[150]: Period('2007', 'A-DEC')

这里,这个Period对象表示的是从2007年1月1日到2007年12月31日之间的整段时间。只需对Period对象加上或减去一个整数即可达到根据其频率进行位移的效果:

In [151]: p   5
Out[151]: Period('2012', 'A-DEC')

In [152]: p - 2
Out[152]: Period('2005', 'A-DEC')

如果两个Period对象拥有相同的频率,则它们的差就是它们之间的单位数量:

In [153]: pd.Period('2014', freq='A-DEC') - p
Out[153]: 7

period_range函数可用于创建规则的时期范围:

In [154]: rng = pd.period_range('2000-01-01', '2000-06-30', freq='M')

In [155]: rng
Out[155]: PeriodIndex(['2000-01', '2000-02', '2000-03', '2000-04', '2000-05', '20
00-06'], dtype='period[M]', freq='M')

PeriodIndex类保存了一组Period,它可以在任何pandas数据结构中被用作轴索引:

In [156]: pd.Series(np.random.randn(6), index=rng)
Out[156]: 
2000-01   -0.514551
2000-02   -0.559782
2000-03   -0.783408
2000-04   -1.797685
2000-05   -0.172670
2000-06    0.680215
Freq: M, dtype: float64

如果你有一个字符串数组,你也可以使用PeriodIndex类:

In [157]: values = ['2001Q3', '2002Q2', '2003Q1']

In [158]: index = pd.PeriodIndex(values, freq='Q-DEC')

In [159]: index
Out[159]: PeriodIndex(['2001Q3', '2002Q2', '2003Q1'], dtype='period[Q-DEC]', freq
='Q-DEC')

时期的频率转换

Period和PeriodIndex对象都可以通过其asfreq方法被转换成别的频率。假设我们有一个年度时期,希望将其转换为当年年初或年末的一个月度时期。该任务非常简单:

In [160]: p = pd.Period('2007', freq='A-DEC')

In [161]: p
Out[161]: Period('2007', 'A-DEC')

In [162]: p.asfreq('M', how='start')
Out[162]: Period('2007-01', 'M')

In [163]: p.asfreq('M', how='end')
Out[163]: Period('2007-12', 'M')

你可以将Period('2007','A-DEC')看做一个被划分为多个月度时期的时间段中的游标。图11-1对此进行了说明。对于一个不以12月结束的财政年度,月度子时期的归属情况就不一样了:

In [164]: p = pd.Period('2007', freq='A-JUN')

In [165]: p
Out[165]: Period('2007', 'A-JUN')

In [166]: p.asfreq('M', 'start')
Out[166]: Period('2006-07', 'M')

In [167]: p.asfreq('M', 'end')
Out[167]: Period('2007-06', 'M')

澳门新浦京娱乐场网站 35

图11-1 Period频率转换示例

在将高频率转换为低频率时,超时期(superperiod)是由子时期(subperiod)所属的位置决定的。例如,在A-JUN频率中,月份“2007年8月”实际上是属于周期“2008年”的:

In [168]: p = pd.Period('Aug-2007', 'M')

In [169]: p.asfreq('A-JUN')
Out[169]: Period('2008', 'A-JUN')

完整的PeriodIndex或TimeSeries的频率转换方式也是如此:

In [170]: rng = pd.period_range('2006', '2009', freq='A-DEC')

In [171]: ts = pd.Series(np.random.randn(len(rng)), index=rng)

In [172]: ts
Out[172]: 
2006    1.607578
2007    0.200381
2008   -0.834068
2009   -0.302988
Freq: A-DEC, dtype: float64

In [173]: ts.asfreq('M', how='start')
Out[173]: 
2006-01    1.607578
2007-01    0.200381
2008-01   -0.834068
2009-01   -0.302988
Freq: M, dtype: float64

这里,根据年度时期的第一个月,每年的时期被取代为每月的时期。如果我们想要每年的最后一个工作日,我们可以使用“B”频率,并指明想要该时期的末尾:

In [174]: ts.asfreq('B', how='end')

Out[174]: 
2006-12-29    1.607578
2007-12-31    0.200381
2008-12-31   -0.834068
2009-12-31   -0.302988
Freq: B, dtype: float64

按季度计算的时期频率

季度型数据在会计、金融等领域中很常见。许多季度型数据都会涉及“财年末”的概念,通常是一年12个月中某月的最后一个日历日或工作日。就这一点来说,时期"2012Q4"根据财年末的不同会有不同的含义。pandas支持12种可能的季度型频率,即Q-JAN到Q-DEC:

In [175]: p = pd.Period('2012Q4', freq='Q-JAN')

In [176]: p
Out[176]: Period('2012Q4', 'Q-JAN')

在以1月结束的财年中,2012Q4是从11月到1月(将其转换为日型频率就明白了)。图11-2对此进行了说明:

In [177]: p.asfreq('D', 'start')
Out[177]: Period('2011-11-01', 'D')

In [178]: p.asfreq('D', 'end')
Out[178]: Period('2012-01-31', 'D')

澳门新浦京娱乐场网站 36

澳门新浦京娱乐场网站时间序列,时间序列分析工具箱。图11.2 不同季度型频率之间的转换

因此,Period之间的算术运算会非常简单。例如,要获取该季度倒数第二个工作日下午4点的时间戳,你可以这样:

In [179]: p4pm = (p.asfreq('B', 'e') - 1).asfreq('T', 's')   16 * 60

In [180]: p4pm
Out[180]: Period('2012-01-30 16:00', 'T')

In [181]: p4pm.to_timestamp()
Out[181]: Timestamp('2012-01-30 16:00:00')

period_range可用于生成季度型范围。季度型范围的算术运算也跟上面是一样的:

In [182]: rng = pd.period_range('2011Q3', '2012Q4', freq='Q-JAN')

In [183]: ts = pd.Series(np.arange(len(rng)), index=rng)

In [184]: ts
Out[184]: 
2011Q3    0
2011Q4    1
2012Q1    2
2012Q2    3
2012Q3    4
2012Q4    5
Freq: Q-JAN, dtype: int64

In [185]: new_rng = (rng.asfreq('B', 'e') - 1).asfreq('T', 's')   16 * 60

In [186]: ts.index = new_rng.to_timestamp()

In [187]: ts
Out[187]:
2010-10-28 16:00:00    0
2011-01-28 16:00:00    1
2011-04-28 16:00:00    2
2011-07-28 16:00:00    3
2011-10-28 16:00:00    4
2012-01-30 16:00:00    5
dtype: int64

将Timestamp转换为Period(及其反向过程)

通过使用to_period方法,可以将由时间戳索引的Series和DataFrame对象转换为以时期索引:

In [188]: rng = pd.date_range('2000-01-01', periods=3, freq='M')

In [189]: ts = pd.Series(np.random.randn(3), index=rng)

In [190]: ts
Out[190]: 
2000-01-31    1.663261
2000-02-29   -0.996206
2000-03-31    1.521760
Freq: M, dtype: float64

In [191]: pts = ts.to_period()

In [192]: pts
Out[192]: 
2000-01    1.663261
2000-02   -0.996206
2000-03    1.521760
Freq: M, dtype: float64

由于时期指的是非重叠时间区间,因此对于给定的频率,一个时间戳只能属于一个时期。新PeriodIndex的频率默认是从时间戳推断而来的,你也可以指定任何别的频率。结果中允许存在重复时期:

In [193]: rng = pd.date_range('1/29/2000', periods=6, freq='D')

In [194]: ts2 = pd.Series(np.random.randn(6), index=rng)

In [195]: ts2
Out[195]: 
2000-01-29    0.244175
2000-01-30    0.423331
2000-01-31   -0.654040
2000-02-01    2.089154
2000-02-02   -0.060220
2000-02-03   -0.167933
Freq: D, dtype: float64

In [196]: ts2.to_period('M')
Out[196]: 
2000-01    0.244175
2000-01    0.423331
2000-01   -0.654040
2000-02    2.089154
2000-02   -0.060220
2000-02   -0.167933
Freq: M, dtype: float64

要转换回时间戳,使用to_timestamp即可:

In [197]: pts = ts2.to_period()

In [198]: pts
Out[198]: 
2000-01-29    0.244175
2000-01-30    0.423331
2000-01-31   -0.654040
2000-02-01    2.089154
2000-02-02   -0.060220
2000-02-03   -0.167933
Freq: D, dtype: float64

In [199]: pts.to_timestamp(how='end')
Out[199]: 
2000-01-29    0.244175
2000-01-30    0.423331
2000-01-31   -0.654040
2000-02-01    2.089154
2000-02-02   -0.060220
2000-02-03   -0.167933
Freq: D, dtype: float64

通过数组创建PeriodIndex

固定频率的数据集通常会将时间信息分开存放在多个列中。例如,在下面这个宏观经济数据集中,年度和季度就分别存放在不同的列中:

In [200]: data = pd.read_csv('examples/macrodata.csv')

In [201]: data.head(5)
Out[201]: 
     year  quarter   realgdp  realcons  realinv  realgovt  realdpi    cpi  
0  1959.0      1.0  2710.349    1707.4  286.898   470.045   1886.9  28.98   
1  1959.0      2.0  2778.801    1733.7  310.859   481.301   1919.7  29.15   
2  1959.0      3.0  2775.488    1751.8  289.226   491.260   1916.4  29.35   
3  1959.0      4.0  2785.204    1753.7  299.356   484.052   1931.3  29.37   
4  1960.0      1.0  2847.699    1770.5  331.722   462.199   1955.5  29.54   
      m1  tbilrate  unemp      pop  infl  realint  
0  139.7      2.82    5.8  177.146  0.00     0.00  
1  141.7      3.08    5.1  177.830  2.34     0.74  
2  140.5      3.82    5.3  178.657  2.74     1.09  
3  140.0      4.33    5.6  179.386  0.27     4.06  
4  139.6      3.50    5.2  180.007  2.31     1.19  

In [202]: data.year
Out[202]: 
0      1959.0
1      1959.0
2      1959.0
3      1959.0
4      1960.0
5      1960.0
6      1960.0
7      1960.0
8      1961.0
9      1961.0
        ...  
193    2007.0
194    2007.0
195    2007.0
196    2008.0
197    2008.0
198    2008.0
199    2008.0
200    2009.0
201    2009.0
202    2009.0
Name: year, Length: 203, dtype: float64

In [203]: data.quarter
Out[203]: 
0      1.0
1      2.0
2      3.0
3      4.0
4      1.0
5      2.0
6      3.0
7      4.0
8      1.0
9      2.0
      ... 
193    2.0
194    3.0
195    4.0
196    1.0
197    2.0
198    3.0
199    4.0
200    1.0
201    2.0
202    3.0
Name: quarter, Length: 203, dtype: float64

通过通过将这些数组以及一个频率传入PeriodIndex,就可以将它们合并成DataFrame的一个索引:

In [204]: index = pd.PeriodIndex(year=data.year, quarter=data.quarter,
   .....:                        freq='Q-DEC')

In [205]: index
Out[205]: 
PeriodIndex(['1959Q1', '1959Q2', '1959Q3', '1959Q4', '1960Q1', '1960Q2',
             '1960Q3', '1960Q4', '1961Q1', '1961Q2',
             ...
             '2007Q2', '2007Q3', '2007Q4', '2008Q1', '2008Q2', '2008Q3',
             '2008Q4', '2009Q1', '2009Q2', '2009Q3'],
            dtype='period[Q-DEC]', length=203, freq='Q-DEC')

In [206]: data.index = index

In [207]: data.infl
Out[207]: 
1959Q1    0.00
1959Q2    2.34
1959Q3    2.74
1959Q4    0.27
1960Q1    2.31
1960Q2    0.14
1960Q3    2.70
1960Q4    1.21
1961Q1   -0.40
1961Q2    1.47
          ... 
2007Q2    2.75
2007Q3    3.45
2007Q4    6.38
2008Q1    2.82
2008Q2    8.53
2008Q3   -3.16
2008Q4   -8.79
2009Q1    0.94
2009Q2    3.37
2009Q3    3.56
Freq: Q-DEC, Name: infl, Length: 203, dtype: float64

11.6 重采样及频率转换

重采样(resampling)指的是将时间序列从一个频率转换到另一个频率的处理过程。将高频率数据聚合到低频率称为降采样(downsampling),而将低频率数据转换到高频率则称为升采样(upsampling)。并不是所有的重采样都能被划分到这两个大类中。例如,将W-WED(每周三)转换为W-FRI既不是降采样也不是升采样。

pandas对象都带有一个resample方法,它是各种频率转换工作的主力函数。resample有一个类似于groupby的API,调用resample可以分组数据,然后会调用一个聚合函数:

In [208]: rng = pd.date_range('2000-01-01', periods=100, freq='D')

In [209]: ts = pd.Series(np.random.randn(len(rng)), index=rng)

In [210]: ts
Out[210]: 
2000-01-01    0.631634
2000-01-02   -1.594313
2000-01-03   -1.519937
2000-01-04    1.108752
2000-01-05    1.255853
2000-01-06   -0.024330
2000-01-07   -2.047939
2000-01-08   -0.272657
2000-01-09   -1.692615
2000-01-10    1.423830
                ...   
2000-03-31   -0.007852
2000-04-01   -1.638806
2000-04-02    1.401227
2000-04-03    1.758539
2000-04-04    0.628932
2000-04-05   -0.423776
2000-04-06    0.789740
2000-04-07    0.937568
2000-04-08   -2.253294
2000-04-09   -1.772919
Freq: D, Length: 100, dtype: float64

In [211]: ts.resample('M').mean()
Out[211]: 
2000-01-31   -0.165893
2000-02-29    0.078606
2000-03-31    0.223811
2000-04-30   -0.063643
Freq: M, dtype: float64

In [212]: ts.resample('M', kind='period').mean()
Out[212]: 
2000-01   -0.165893
2000-02    0.078606
2000-03    0.223811
2000-04   -0.063643
Freq: M, dtype: float64

resample是一个灵活高效的方法,可用于处理非常大的时间序列。我将通过一系列的示例说明其用法。表11-5总结它的一些选项。

表11-5 resample方法的参数

澳门新浦京娱乐场网站 37

降采样

将数据聚合到规律的低频率是一件非常普通的时间序列处理任务。待聚合的数据不必拥有固定的频率,期望的频率会自动定义聚合的面元边界,这些面元用于将时间序列拆分为多个片段。例如,要转换到月度频率('M'或'BM'),数据需要被划分到多个单月时间段中。各时间段都是半开放的。一个数据点只能属于一个时间段,所有时间段的并集必须能组成整个时间帧。在用resample对数据进行降采样时,需要考虑两样东西:

  • 各区间哪边是闭合的。
  • 如何标记各个聚合面元,用区间的开头还是末尾。

为了说明,我们来看一些“1分钟”数据:

In [213]: rng = pd.date_range('2000-01-01', periods=12, freq='T')

In [214]: ts = pd.Series(np.arange(12), index=rng)

In [215]: ts
Out[215]: 
2000-01-01 00:00:00     0
2000-01-01 00:01:00     1
2000-01-01 00:02:00     2
2000-01-01 00:03:00     3
2000-01-01 00:04:00     4
2000-01-01 00:05:00     5
2000-01-01 00:06:00     6
2000-01-01 00:07:00     7
2000-01-01 00:08:00     8
2000-01-01 00:09:00     9
2000-01-01 00:10:00    10
2000-01-01 00:11:00    11
Freq: T, dtype: int64

假设你想要通过求和的方式将这些数据聚合到“5分钟”块中:

In [216]: ts.resample('5min', closed='right').sum()
Out[216]: 
1999-12-31 23:55:00     0
2000-01-01 00:00:00    15
2000-01-01 00:05:00    40
2000-01-01 00:10:00    11
Freq: 5T, dtype: int64

传入的频率将会以“5分钟”的增量定义面元边界。默认情况下,面元的右边界是包含的,因此00:00到00:05的区间中是包含00:05的。传入closed='left'会让区间以左边界闭合:

In [217]: ts.resample('5min', closed='right').sum()
Out[217]: 
1999-12-31 23:55:00     0
2000-01-01 00:00:00    15
2000-01-01 00:05:00    40
2000-01-01 00:10:00    11
Freq: 5T, dtype: int64

如你所见,最终的时间序列是以各面元右边界的时间戳进行标记的。传入label='right'即可用面元的邮编界对其进行标记:

In [218]: ts.resample('5min', closed='right', label='right').sum()
Out[218]: 
2000-01-01 00:00:00     0
2000-01-01 00:05:00    15
2000-01-01 00:10:00    40
2000-01-01 00:15:00    11
Freq: 5T, dtype: int64

图11-3说明了“1分钟”数据被转换为“5分钟”数据的处理过程。

澳门新浦京娱乐场网站 38

图11-3 各种closed、label约定的“5分钟”重采样演示

最后,你可能希望对结果索引做一些位移,比如从右边界减去一秒以便更容易明白该时间戳到底表示的是哪个区间。只需通过loffset设置一个字符串或日期偏移量即可实现这个目的:

In [219]: ts.resample('5min', closed='right',
   .....:             label='right', loffset='-1s').sum()
Out[219]: 
1999-12-31 23:59:59     0
2000-01-01 00:04:59    15
In [219]: ts.resample('5min', closed='right',
   .....:             label='right', loffset='-1s').sum()
Out[219]: 
1999-12-31 23:59:59     0
2000-01-01 00:04:59    15

此外,也可以通过调用结果对象的shift方法来实现该目的,这样就不需要设置loffset了。

OHLC重采样

金融领域中有一种无所不在的时间序列聚合方式,即计算各面元的四个值:第一个值(open,开盘)、最后一个值(close,收盘)、最大值(high,最高)以及最小值(low,最低)。传入how='ohlc'即可得到一个含有这四种聚合值的DataFrame。整个过程很高效,只需一次扫描即可计算出结果:

In [220]: ts.resample('5min').ohlc()
Out[220]: 
                     open  high  low  close
2000-01-01 00:00:00     0     4    0      4
2000-01-01 00:05:00     5     9    5      9
2000-01-01 00:10:00    10    11   10     11

升采样和插值

在将数据从低频率转换到高频率时,就不需要聚合了。我们来看一个带有一些周型数据的DataFrame:

In [221]: frame = pd.DataFrame(np.random.randn(2, 4),
   .....:                      index=pd.date_range('1/1/2000', periods=2,
   .....:                                          freq='W-WED'),
   .....:                      columns=['Colorado', 'Texas', 'New York', 'Ohio'])

In [222]: frame
Out[222]: 
            Colorado     Texas  New York      Ohio
2000-01-05 -0.896431  0.677263  0.036503  0.087102
2000-01-12 -0.046662  0.927238  0.482284 -0.867130

当你对这个数据进行聚合,每组只有一个值,这样就会引入缺失值。我们使用asfreq方法转换成高频,不经过聚合:

In [223]: df_daily = frame.resample('D').asfreq()

In [224]: df_daily
Out[224]: 
            Colorado     Texas  New York      Ohio
2000-01-05 -0.896431  0.677263  0.036503  0.087102
2000-01-06       NaN       NaN       NaN       NaN
2000-01-07       NaN       NaN       NaN       NaN
2000-01-08       NaN       NaN       NaN       NaN
2000-01-09       NaN       NaN       NaN       NaN
2000-01-10       NaN       NaN       NaN       NaN
2000-01-11       NaN       NaN       NaN       NaN
2000-01-12 -0.046662  0.927238  0.482284 -0.867130

假设你想要用前面的周型值填充“非星期三”。resampling的填充和插值方式跟fillna和reindex的一样:

In [225]: frame.resample('D').ffill()
Out[225]: 
            Colorado     Texas  New York      Ohio
2000-01-05 -0.896431  0.677263  0.036503  0.087102
2000-01-06 -0.896431  0.677263  0.036503  0.087102
2000-01-07 -0.896431  0.677263  0.036503  0.087102
2000-01-08 -0.896431  0.677263  0.036503  0.087102
2000-01-09 -0.896431  0.677263  0.036503  0.087102
2000-01-10 -0.896431  0.677263  0.036503  0.087102
2000-01-11 -0.896431  0.677263  0.036503  0.087102
2000-01-12 -0.046662  0.927238  0.482284 -0.867130

同样,这里也可以只填充指定的时期数(目的是限制前面的观测值的持续使用距离):

In [226]: frame.resample('D').ffill(limit=2)
Out[226]:
            Colorado     Texas  New York      Ohio
2000-01-05 -0.896431  0.677263  0.036503  0.087102
2000-01-06 -0.896431  0.677263  0.036503  0.087102
2000-01-07 -0.896431  0.677263  0.036503  0.087102
2000-01-08       NaN       NaN       NaN       NaN
2000-01-09       NaN       NaN       NaN       NaN
2000-01-10       NaN       NaN       NaN       NaN
2000-01-11       NaN       NaN       NaN       NaN
2000-01-12 -0.046662  0.927238  0.482284 -0.867130

注意,新的日期索引完全没必要跟旧的重叠:

In [227]: frame.resample('W-THU').ffill()
Out[227]: 
            Colorado     Texas  New York      Ohio
2000-01-06 -0.896431  0.677263  0.036503  0.087102
2000-01-13 -0.046662  0.927238  0.482284 -0.867130

通过时期进行重采样

对那些使用时期索引的数据进行重采样与时间戳很像:

In [228]: frame = pd.DataFrame(np.random.randn(24, 4),
   .....:                      index=pd.period_range('1-2000', '12-2001',
   .....:                                            freq='M'),
   .....:                      columns=['Colorado', 'Texas', 'New York', 'Ohio'])

In [229]: frame[:5]
Out[229]: 
         Colorado     Texas  New York      Ohio
2000-01  0.493841 -0.155434  1.397286  1.507055
2000-02 -1.179442  0.443171  1.395676 -0.529658
2000-03  0.787358  0.248845  0.743239  1.267746
2000-04  1.302395 -0.272154 -0.051532 -0.467740
2000-05 -1.040816  0.426419  0.312945 -1.115689

In [230]: annual_frame = frame.resample('A-DEC').mean()

In [231]: annual_frame
Out[231]: 
      Colorado     Texas  New York      Ohio
2000  0.556703  0.016631  0.111873 -0.027445
2001  0.046303  0.163344  0.251503 -0.157276

升采样要稍微麻烦一些,因为你必须决定在新频率中各区间的哪端用于放置原来的值,就像asfreq方法那样。convention参数默认为'start',也可设置为'end':

# Q-DEC: Quarterly, year ending in December
In [232]: annual_frame.resample('Q-DEC').ffill()
Out[232]: 
        Colorado     Texas  New York      Ohio
2000Q1  0.556703  0.016631  0.111873 -0.027445
2000Q2  0.556703  0.016631  0.111873 -0.027445
2000Q3  0.556703  0.016631  0.111873 -0.027445
2000Q4  0.556703  0.016631  0.111873 -0.027445
2001Q1  0.046303  0.163344  0.251503 -0.157276
2001Q2  0.046303  0.163344  0.251503 -0.157276
2001Q3  0.046303  0.163344  0.251503 -0.157276
2001Q4  0.046303  0.163344  0.251503 -0.157276

In [233]: annual_frame.resample('Q-DEC', convention='end').ffill()
Out[233]: 
        Colorado     Texas  New York      Ohio
2000Q4  0.556703  0.016631  0.111873 -0.027445
2001Q1  0.556703  0.016631  0.111873 -0.027445
2001Q2  0.556703  0.016631  0.111873 -0.027445
2001Q3  0.556703  0.016631  0.111873 -0.027445
2001Q4  0.046303  0.163344  0.251503 -0.157276

由于时期指的是时间区间,所以升采样和降采样的规则就比较严格:

  • 在降采样中,目标频率必须是源频率的子时期(subperiod)。
  • 在升采样中,目标频率必须是源频率的超时期(superperiod)。

如果不满足这些条件,就会引发异常。这主要影响的是按季、年、周计算的频率。例如,由Q-MAR定义的时间区间只能升采样为A-MAR、A-JUN、A-SEP、A-DEC等:

In [234]: annual_frame.resample('Q-MAR').ffill()
Out[234]: 
        Colorado     Texas  New York      Ohio
2000Q4  0.556703  0.016631  0.111873 -0.027445
2001Q1  0.556703  0.016631  0.111873 -0.027445
2001Q2  0.556703  0.016631  0.111873 -0.027445
2001Q3  0.556703  0.016631  0.111873 -0.027445
2001Q4  0.046303  0.163344  0.251503 -0.157276
2002Q1  0.046303  0.163344  0.251503 -0.157276
2002Q2  0.046303  0.163344  0.251503 -0.157276
2002Q3  0.046303  0.163344  0.251503 -0.157276

11.7 移动窗口函数

在移动窗口(可以带有指数衰减权数)上计算的各种统计函数也是一类常见于时间序列的数组变换。这样可以圆滑噪音数据或断裂数据。我将它们称为移动窗口函数(moving window function),其中还包括那些窗口不定长的函数(如指数加权移动平均)。跟其他统计函数一样,移动窗口函数也会自动排除缺失值。

开始之前,我们加载一些时间序列数据,将其重采样为工作日频率:

In [235]: close_px_all = pd.read_csv('examples/stock_px_2.csv',
   .....:                            parse_dates=True, index_col=0)

In [236]: close_px = close_px_all[['AAPL', 'MSFT', 'XOM']]

In [237]: close_px = close_px.resample('B').ffill()

现在引入rolling运算符,它与resample和groupby很像。可以在TimeSeries或DataFrame以及一个window(表示期数,见图11-4)上调用它:

In [238]: close_px.AAPL.plot()
Out[238]: <matplotlib.axes._subplots.AxesSubplot at 0x7f2f2570cf98>

In [239]: close_px.AAPL.rolling(250).mean().plot()

澳门新浦京娱乐场网站 39

图11-4 苹果公司股价的250日均线

表达式rolling(250)与groupby很像,但不是对其进行分组,而是创建一个按照250天分组的滑动窗口对象。然后,我们就得到了苹果公司股价的250天的移动窗口。

默认情况下,rolling函数需要窗口中所有的值为非NA值。可以修改该行为以解决缺失数据的问题。其实,在时间序列开始处尚不足窗口期的那些数据就是个特例(见图11-5):

In [241]: appl_std250 = close_px.AAPL.rolling(250, min_periods=10).std()

In [242]: appl_std250[5:12]
Out[242]: 
2003-01-09         NaN
2003-01-10         NaN
2003-01-13         NaN
2003-01-14         NaN
2003-01-15    0.077496
2003-01-16    0.074760
2003-01-17    0.112368
Freq: B, Name: AAPL, dtype: float64

In [243]: appl_std250.plot()

澳门新浦京娱乐场网站 40

图11-5 苹果公司250日每日回报标准差

要计算扩展窗口平均(expanding window mean),可以使用expanding而不是rolling。“扩展”意味着,从时间序列的起始处开始窗口,增加窗口直到它超过所有的序列。apple_std250时间序列的扩展窗口平均如下所示:

In [244]: expanding_mean = appl_std250.expanding().mean()

对DataFrame调用rolling_mean(以及与之类似的函数)会将转换应用到所有的列上(见图11-6):

In [246]: close_px.rolling(60).mean().plot(logy=True)

澳门新浦京娱乐场网站 41

图11-6 各股价60日均线(对数Y轴)

rolling函数也可以接受一个指定固定大小时间补偿字符串,而不是一组时期。这样可以方便处理不规律的时间序列。这些字符串也可以传递给resample。例如,我们可以计算20天的滚动均值,如下所示:

In [247]: close_px.rolling('20D').mean()
Out[247]:
                  AAPL       MSFT        XOM
2003-01-02    7.400000  21.110000  29.220000
2003-01-03    7.425000  21.125000  29.230000
2003-01-06    7.433333  21.256667  29.473333
2003-01-07    7.432500  21.425000  29.342500
2003-01-08    7.402000  21.402000  29.240000
2003-01-09    7.391667  21.490000  29.273333
2003-01-10    7.387143  21.558571  29.238571
2003-01-13    7.378750  21.633750  29.197500
2003-01-14    7.370000  21.717778  29.194444
2003-01-15    7.355000  21.757000  29.152000
...                ...        ...        ...
2011-10-03  398.002143  25.890714  72.413571
2011-10-04  396.802143  25.807857  72.427143
2011-10-05  395.751429  25.729286  72.422857
2011-10-06  394.099286  25.673571  72.375714
2011-10-07  392.479333  25.712000  72.454667
2011-10-10  389.351429  25.602143  72.527857
2011-10-11  388.505000  25.674286  72.835000
2011-10-12  388.531429  25.810000  73.400714
2011-10-13  388.826429  25.961429  73.905000
2011-10-14  391.038000  26.048667  74.185333
[2292 rows x 3 columns]

指数加权函数

另一种使用固定大小窗口及相等权数观测值的办法是,定义一个衰减因子(decay factor)常量,以便使近期的观测值拥有更大的权数。衰减因子的定义方式有很多,比较流行的是使用时间间隔(span),它可以使结果兼容于窗口大小等于时间间隔的简单移动窗口(simple moving window)函数。

由于指数加权统计会赋予近期的观测值更大的权数,因此相对于等权统计,它能“适应”更快的变化。

除了rolling和expanding,pandas还有ewm运算符。下面这个例子对比了苹果公司股价的30日移动平均和span=30的指数加权移动平均(如图11-7所示):

In [249]: aapl_px = close_px.AAPL['2006':'2007']

In [250]: ma60 = aapl_px.rolling(30, min_periods=20).mean()

In [251]: ewma60 = aapl_px.ewm(span=30).mean()

In [252]: ma60.plot(style='k--', label='Simple MA')
Out[252]: <matplotlib.axes._subplots.AxesSubplot at 0x7f2f252161d0>

In [253]: ewma60.plot(style='k-', label='EW MA')
Out[253]: <matplotlib.axes._subplots.AxesSubplot at 0x7f2f252161d0>

In [254]: plt.legend()

澳门新浦京娱乐场网站 42

图11-7 简单移动平均与指数加权移动平均

二元移动窗口函数

有些统计运算(如相关系数和协方差)需要在两个时间序列上执行。例如,金融分析师常常对某只股票对某个参考指数(如标准普尔500指数)的相关系数感兴趣。要进行说明,我们先计算我们感兴趣的时间序列的百分数变化:

In [256]: spx_px = close_px_all['SPX']

In [257]: spx_rets = spx_px.pct_change()

In [258]: returns = close_px.pct_change()

调用rolling之后,corr聚合函数开始计算与spx_rets滚动相关系数(结果见图11-8):

In [259]: corr = returns.AAPL.rolling(125, min_periods=100).corr(spx_rets)

In [260]: corr.plot()

澳门新浦京娱乐场网站 43

图11-8 AAPL 6个月的回报与标准普尔500指数的相关系数

假设你想要一次性计算多只股票与标准普尔500指数的相关系数。虽然编写一个循环并新建一个DataFrame不是什么难事,但比较啰嗦。其实,只需传入一个TimeSeries和一个DataFrame,rolling_corr就会自动计算TimeSeries(本例中就是spx_rets)与DataFrame各列的相关系数。结果如图11-9所示:

In [262]: corr = returns.rolling(125, min_periods=100).corr(spx_rets)

In [263]: corr.plot()

澳门新浦京娱乐场网站 44

图11-9 3只股票6个月的回报与标准普尔500指数的相关系数

用户定义的移动窗口函数

rolling_apply函数使你能够在移动窗口上应用自己设计的数组函数。唯一要求的就是:该函数要能从数组的各个片段中产生单个值(即约简)。比如说,当我们用rolling(...).quantile(q)计算样本分位数时,可能对样本中特定值的百分等级感兴趣。scipy.stats.percentileofscore函数就能达到这个目的(结果见图11-10):

In [265]: from scipy.stats import percentileofscore

In [266]: score_at_2percent = lambda x: percentileofscore(x, 0.02)

In [267]: result = returns.AAPL.rolling(250).apply(score_at_2percent)

In [268]: result.plot()

澳门新浦京娱乐场网站 45

图11-10 AAPL 2%回报率的百分等级(一年窗口期)

如果你没安装SciPy,可以使用conda或pip安装。

11.8 总结

与前面章节接触的数据相比,时间序列数据要求不同类型的分析和数据转换工具。

在接下来的章节中,我们将学习一些高级的pandas方法和如何开始使用建模库statsmodels和scikit-learn。


第1章 准备工作
第2章 Python语法基础,IPython和Jupyter
第3章 Python的数据结构、函数和文件
第4章 NumPy基础:数组和矢量计算
第5章 pandas入门
第6章 数据加载、存储与文件格式
第7章 数据清洗和准备
第8章 数据规整:聚合、合并和重塑
第9章 绘图和可视化
第10章 数据聚合与分组运算
第11章 时间序列
第12章 pandas高级应用
第13章 Python建模库介绍
第14章 数据分析案例
附录A NumPy高级应用
附录B 更多关于IPython的内容(完)


本文由澳门新浦京娱乐场网站发布于www.146.net,转载请注明出处:澳门新浦京娱乐场网站时间序列,时间序列分析