日历中的夏天

日历中的夏天 | 统计之都

看着有节，摸着无节，打一生活用品

不知不觉，夏日已慢慢临近。姑娘们飞扬的裙角，小贩叫卖的西瓜，蚊蝇嗡嗡的声音，以及翻过的一页日历，都提醒着你夏天快来了。夏季有着不同的定义：根据中国人的日历，我们所俗称的夏季从立夏开始，到立秋结束；但在气候学上，夏季是指连续五天平均温度超过 22 摄氏度即算作夏季的开始，若连续五天平均温度低于 22 度则算作入秋；而天文学上的夏季一般是指六、七、八这三个月。

那么哪一种夏季的定义更合适一些呢？还是用数据可视化来说话吧。这项任务基本上有两个步骤：一是获取某城市的 2011 年日平均温度数据，二是根据数据绘制日历热图（Calendar-Heatmap）。

本文所采用的数据源是 Wunderground 提供的 API。该 API 所提供的数据极为丰富，除了历史温度数据之外，还有湿度、风向等大量信息可供利用，所以它也被数据源手册所推荐。为了使用这个 API 我们先申请一个免费帐户，然后利用 R 语言中的 RJSONIO 包来提取每天的平均温度。使用这个 API 要注意的是，其免费帐户限制了每分钟 10 次调用，超出会中断连接。不知道别人怎么弄的，本人的笨办法是在程序中增加了暂停。这样获得了 365 个日平均温度。

日历热图是一种有趣的工具，它可以在日历表中显示时间序列数据的变化。在 2009 年的 The Data Expo 中，获奖团队就是利用 SAS 来生成日历热图。多才多艺的 R 语言当然也可以做到。在《R Graphs Cookbook》这本书中就提到了绘制日历热图的方法，第一种方法是 Paul Bleicher 所写的一个函数，它是基于 grid,lattice,chron 这三个扩展包来编写的。第二种方法是使用 openair 扩展包中的 calendarPlot 函数。生成的图形就象下面这个样子。

看起来不错，但是我们还没完。我们希望挑选出平均温度在 22 度以上的日子，突出显示出来。所以我们采用第三种方法，用 ggplot2 包来绘制日历热图，图形显示如下（参考了 MarginTale 的这篇文章）。

上图数日子是竖着来数的，横轴表示每月的第几周，纵轴表示星期几。灰色部分表示当天平均温度在 22 度以下，有色彩的区块表示在 22 度以上。颜色越偏黄则表示温度越高。在 2011 年，立夏的时间是 5 月 6 日，立秋是 8 月 8 日，但可以观察到立秋之后仍有很多日子的平均温度在 22 度以上。这就是我们所俗称的 “秋老虎”。如果按照气候学的定义，四月末就有五天以上连续的高温天气，照这样看夏天应该在四月末就开始了，一直延续到十月初结束。而天文学上的夏季则是六、七、八三个月，看到这三个月基本上全是 22 度以上，而且高温天气集中在七八两个月，这也正是学校放暑假的时间段。这样看来，似乎天文学的夏季定义是比较符合我们人体的感觉的。其它的要么偏短，要么偏长。

写到这里，想起了梁静茹的一首歌：

宁静的夏天

天空中繁星点点

心里头有些思念

思念着你的脸

我可以假装看不见

也可以偷偷的想念

直到让我摸到你那温暖的脸

（最后要说的是，本人并非气象专家，本文也没有考虑到湿度对体感温度的影响，或是其它因素。主要还是向各位介绍 R 语言中日历热图的绘制以及数据的获取。）

代码如下：

# 加载所需扩展包
library(RCurl)
library(RJSONIO)
require(quantmod)
library(ggplot2)

# 提取武汉市2011年一年的历史数据
date<- seq.Date(from=as.Date('2011-01-01'),
  to=as.Date('2011-12-31'), by='1 day')
date.range<- as.character(format(date,"%Y%m%d"))
n <- length(date.range)
temp <- humi <- rep(0,n)
for (i in 1:n) {
  # 你要用自己申请的API key来代替程序中的yourkey
  url <- 'http://api.wunderground.com/api/yourkey/'
  finalurl <- paste(url,'history_',date.range[i],
    '/q/wuhan.json',sep='')
  web <- getURL(finalurl)
  raw <-fromJSON(web)
  temp[i] <- raw$history$dailysummary[[1]]$meantempm
  humi[i] <- raw$history$dailysummary[[1]]$humidity
  # 在循环内增加一个7秒的暂停，避免连接断开。
  Sys.sleep(7)
}
# 将获得的数据整合为数据框，并将温度和湿度转为数值格式
dataset <- data.frame(temp,humi,date,stringsAsFactors=F)
dataset$temp <- as.numeric(dataset$temp)
dataset$humi <- as.numeric(dataset$humi)

# 用openair包绘制日历热图
install.packages('openair')
library(openair)
calendarPlot(dataset,pollutant='temp',year=2011)

# 用ggplot2包绘制日历热图
# 复制一个新的数据框
dat <- dataset
# 先取得月份，再转为因子格式
dat$month<-as.numeric(as.POSIXlt(dat$date)$mon+1)
dat$monthf<-factor(dat$month,levels=as.character(1:12),
  labels=c("Jan","Feb","Mar","Apr","May","Jun","Jul",
  "Aug","Sep","Oct","Nov","Dec"),ordered=TRUE)
# 得到每周的星期，也转为因子格式
dat$weekday = as.POSIXlt(dat$date)$wday
dat$weekdayf<-factor(dat$weekday,levels=rev(0:6),
  labels=rev(c("Sun","Mon","Tue","Wed","Thu","Fri","Sat")),ordered=TRUE)
# 先得到全年的周序号，然后得到每个月的周序号
dat$week <- as.numeric(format(dat$date,"%W"))
dat<-ddply(dat,.(monthf),transform,monthweek=1+week-min(week))
# 绘图
P <- ggplot(dat, aes(monthweek, weekdayf, fill = temp)) +
  geom_tile(colour='white') +
  facet_wrap(~monthf ,nrow=3) +
  scale_fill_gradient(space="Lab",limits=c(22, max(dat$value)),
    low="red", high="yellow") +
  opts(title = "武汉市2011年气温日历热图") +
  xlab("Week of Month") + ylab("")
P

一个喜欢折腾数据的人，与李舰合著了《数据科学中的 R 语言》，现就职于蚂蚁金服。