過去の外来患者数から今後の患者数を予測してみる

  1. 電子カルテデータベースから読み込んだ、外来患者数ファイルを取り込み。
  2. xtsパッケージで時系列データの取り扱い。
  3. データの視覚化。プロットしてみる。
  4. 時系列解析。

予め、読み込み済みのデータを使います。

library(dplyr)
library(ggplot2)
ds %>% tail(n=3)
##          date sinkann saisin zenkannjya jitukannjya
## 51 2015-02-01      97    948       1045         853
## 52 2015-03-01     104   1065       1169         910
## 53 2015-04-01     101   1028       1129         909

元データは、データフレーム型です。 teramonagi先生が、xtsが便利と申されていたので…
ここを参考に
http://www.slideshare.net/teramonagi/tokyo-r15-20110702

library(zoo)
library(xts)
ds1 <- as.xts(read.zoo(ds))
ds1 %>% class
## [1] "xts" "zoo"
ds1 %>% tail(n=3)
##            sinkann saisin zenkannjya jitukannjya
## 2015-02-01      97    948       1045         853
## 2015-03-01     104   1065       1169         910
## 2015-04-01     101   1028       1129         909

xts変換後は、rownamesが日付になっています。

ds1["2014-04-01::2014-06-01"]
##            sinkann saisin zenkannjya jitukannjya
## 2014-04-01     109    780        889         703
## 2014-05-01      90    896        986         776
## 2014-06-01     114    862        976         766

とか、

ds1["2014-04-01::"]
ds1["::2014-04-01"]

とか、日付を指定できます。

欠損値補間とか、年ごとに集計とかいろいろできるそうです。　　

次はプロットしてみましょう。

plot.zoo(ds1, col = c("red", "blue","green","yellow"), plot.type = "single")

ここから時系列解析。

元データの月毎の実患者数をts型に変更。

ds$jitukannjya
##  [1]  183  187  211  206  195  238  251  235  235  325  342  439  409  475
## [15]  414  482  389  397  433  459  432  409  537  669  723  650  583  611
## [29]  621  672  613  589  574  624  738  873  872  845  785  870  703  776
## [43]  766  765  692  746  878  901 1171  861  853  910  909
test<-ts(ds$jitukannjya,start=c(2010,12),frequency=12)
test
##       Jan  Feb  Mar  Apr  May  Jun  Jul  Aug  Sep  Oct  Nov  Dec
## 2010                                                         183
## 2011  187  211  206  195  238  251  235  235  325  342  439  409
## 2012  475  414  482  389  397  433  459  432  409  537  669  723
## 2013  650  583  611  621  672  613  589  574  624  738  873  872
## 2014  845  785  870  703  776  766  765  692  746  878  901 1171
## 2015  861  853  910  909
library(forecast)
seasonplot(test,col=rainbow(12),season.labels = T, year.labels = T,year.labels.left = T,cex=0.8)

まだ開院して５年弱なので、増患傾向かな。

acf(test)

plot(decompose(test))

夏に減患する季節変動は、どこの施設も同じかな…

plot(forecast(test, level = c(95), h = 20))

forecast(test, level = c(95), h = 12)
##          Point Forecast    Lo 95     Hi 95
## May 2015       886.7542 763.3447 1010.1637
## Jun 2015       874.8864 748.0942 1001.6787
## Jul 2015       831.3220 705.0934  957.5506
## Aug 2015       748.7794 629.0570  868.5018
## Sep 2015       844.6148 701.8590  987.3705
## Oct 2015       964.7878 791.9377 1137.6380
## Nov 2015      1086.4250 879.7493 1293.1008
## Dec 2015      1167.6695 931.6017 1403.7373
## Jan 2016      1097.7165 861.8277 1333.6053
## Feb 2016      1018.1178 785.6546 1250.5810
## Mar 2016      1047.4981 793.5695 1301.4266
## Apr 2016       945.1933 702.1872 1188.1995

今年の冬は、1000人超が続きそうなので、看護師さん増やさないと… 　　

時系列解析のお話は、　　　　

ここのサイトがわかりやすいので、参考に。

http://logics-of-blue.com/時系列解析_理論編/

もう少し勉強しないと…
わからないことが多すぎる。
でも、いろいろ使えそうです。