yzhlinscau的个人博客分享 http://blog.sciencenet.cn/u/yzhlinscau

博文

利用R语言简单的时间序列模型预测韩国、日本、伊朗和意大利的疫情趋势

已有 4872 次阅读 2020-3-2 15:09 |个人分类:R|系统分类:科研笔记

问题:韩国、日本、伊朗和意大利的新型冠状病毒肺炎疫情是除了我国之外,目前疫情最严重的4个国家,那么之后疫情会有怎样的变化趋势?
方法:采用R语言forecast包的auto.arima函数进行预测。

1 数据采集
利用腾讯微信提供的数据。仅收集2月28日至今的数据,具体如下:

> df1
# A tibble: 13 x 5
   date                no.kor no.Jpn no.Irn no.Ita
   <dttm>               <dbl>  <dbl>  <dbl>  <dbl>
 1 2020-02-18 00:00:00     31     74      0      3
 2 2020-02-19 00:00:00     51     84      2      3
 3 2020-02-20 00:00:00    104     93      5      3
 4 2020-02-21 00:00:00    208    103     18      6
 5 2020-02-22 00:00:00    433    122     28     20
 6 2020-02-23 00:00:00    602    136     43    132
 7 2020-02-24 00:00:00    833    159     61    230
 8 2020-02-25 00:00:00    977    162     95    283
 9 2020-02-26 00:00:00   1261    178    139    374
10 2020-02-27 00:00:00   1766    207    245    528
11 2020-02-28 00:00:00   2337    233    388    655
12 2020-02-29 00:00:00   3150    241    593    888
13 2020-03-01 00:00:00   3736    256    978   1128

2 批量分析函数的编写
编写下述代码以批量分析,同时预测未来的10天数据:

library(forecast)
#select proper arima model
covpred <- function(dat) {
  yts<-ts(dat$y)
  mod.1 <- auto.arima(yts)
  tt<- forecast(mod.1,10)
  pdf <- data.frame(date = seq(as.Date("2020-03-02"), as.Date("2020-03-11"), "days"), 

                 y = tt$mean[1:10])
 
  return(pdf)
}

3 运行程序
为利用plyr包的批量分析,需要对原始数据做数据重构,过程如下:
df1a<-reshape2::melt(df1,id=1)
head(df1a)
names(df1a)[2:3]<-c('nation','y')

# 批量运行
res <- plyr::dlply(df1a,"nation",covpred)

#预测结果转换
res1<-do.call(cbind,res)
res1a<-res1[,c(1:2,4,6,8)]
names(res1a)<-names(df1)
# 与原始数据合并
df1b<-rbind(df1,res1a)

4 结果查看
批量运行结果如下:
> names(res)
[1] "no.kor" "no.Jpn" "no.Irn" "no.Ita"
> str(res)
List of 4
 $ no.kor:'data.frame': 10 obs. of  2 variables:
  ..$ date: Date[1:10], format: "2020-03-02" "2020-03-03" ...
  ..$ y   : num [1:10] 4322 4908 5494 6080 6666 ...
 $ no.Jpn:'data.frame': 10 obs. of  2 variables:
  ..$ date: Date[1:10], format: "2020-03-02" "2020-03-03" ...
  ..$ y   : num [1:10] 271 286 302 317 332 ...
 $ no.Irn:'data.frame': 10 obs. of  2 variables:
  ..$ date: Date[1:10], format: "2020-03-02" "2020-03-03" ...
  ..$ y   : num [1:10] 1363 1748 2133 2518 2903 ...
 $ no.Ita:'data.frame': 10 obs. of  2 variables:
  ..$ date: Date[1:10], format: "2020-03-02" "2020-03-03" ...
  ..$ y   : num [1:10] 1368 1608 1848 2088 2328 ...
 - attr(*, "split_type")= chr "data.frame"
 - attr(*, "split_labels")='data.frame': 4 obs. of  1 variable:
  ..$ nation: Factor w/ 4 levels "no.kor","no.Jpn",..: 1 2 3 4

5 图形绘制
ggplot2绘图代码如下:

datebreaks <- seq(as.Date("2020-02-18"),as.Date("2020-03-11"),by="3 day")
fnames <- c(
  `no.kor` = "韩国",
  `no.Jpn` = "日本",
  `no.Irn` = "伊朗",
  `no.Ita` = "意大利"
)


ggplot(df1c, aes(as.Date(date), value)) + geom_line() +geom_point()+
  #scale_x_date(format = "%b%Y") + xlab("") + ylab("no.kor")
 geom_text(aes(label=ifelse(as.Date(date)>as.Date('2020-03-01'),as.character(value),'')),
           col='red',hjust=1.15,vjust=.25)+
 geom_text(aes(label=ifelse(as.Date(date)<=as.Date('2020-03-01'),as.character(value),'')),
           col='blue',size=3,vjust=1.5)+
  scale_x_date(breaks=datebreaks)+xlab('日期')+ylab('累计病例')+
  facet_wrap(~variable, scales="free", labeller = as_labeller(fnames))+
  theme(axis.text.x = element_text(angle=45,hjust=1))

生成的图形如下:


从图形来看,韩国和伊朗形势比较严峻,意大利也不乐观,日本总体上比较稳定,爆发的可能性不大。

说明:本文所用的方法属于简单的模型,不考虑数据隐瞒性和人为措施的干预性,因此,模型的可靠性未必高。


————————————————
版权声明:本文为CSDN博主「yzhlinscau」的原创文章,遵循 CC 4.0 BY-SA 版权协议,转载请附上原文出处链接及本声明。
原文链接:https://blog.csdn.net/yzhlinscau/article/details/104611698




https://blog.sciencenet.cn/blog-1114360-1221395.html

上一篇:批量求解基于电导率的植物半致死温度
下一篇:用R分析COVID-19流行病学[译文]
收藏 IP: 219.137.67.*| 热度|

0

该博文允许注册用户评论 请点击登录 评论 (0 个评论)

数据加载中...
扫一扫,分享此博文

Archiver|手机版|科学网 ( 京ICP备07017567号-12 )

GMT+8, 2024-3-28 21:41

Powered by ScienceNet.cn

Copyright © 2007- 中国科学报社

返回顶部