大学入試センター令和７年度試験の問題作成の方向性，試作問題等のなかで、情報Ⅰデータサイエンス分野の参考問題が掲載されています。

  • データサイエンスの有用さが分かる問題であること
  • 興味と知的好奇心をそそる題材であること
  • 将来のデータサイエンスの学びにつながるものであること
  • 日常生活の中でのデータ分析に誘うものであること

これらの要請にかなう問題として作られ、データは公開されている実データが使われています。
日頃からRやPythonなどを利用したデータ分析に親しんでおくことが求められています。

データは、

　一般社団法人日本冷凍空調工業会
　一般社団法人日本アイスクリーム協会
　気象庁 （金沢市のデータ）
からダウンロードまたはコピーし、エクセルで整形しました。

第４問

分析に必要なパッケージを読み込みます。

library(TSstudio)
library(MLmetrics)
## 
##  次のパッケージを付け加えます: 'MLmetrics'
##  以下のオブジェクトは 'package:base' からマスクされています: 
## 
##      Recall
library(forecast)
## Registered S3 method overwritten by 'quantmod':
##   method            from
##   as.zoo.data.frame zoo
library(TTR)
library(quantmod)
##  要求されたパッケージ xts をロード中です
##  要求されたパッケージ zoo をロード中です
## 
##  次のパッケージを付け加えます: 'zoo'
##  以下のオブジェクトは 'package:base' からマスクされています: 
## 
##      as.Date, as.Date.numeric
## 
##  次のパッケージを付け加えます: 'xts'
##  以下のオブジェクトは 'package:dplyr' からマスクされています: 
## 
##      first, last

データの概要

エアコンの売上台数とアイスクリームの売上個数のデータを読み込みます。
気温(kion)、湿度(situdo)は問5で使用します。

x <- read.csv("data/ice_aircon.csv")
str(x)
## 'data.frame':    60 obs. of  5 variables:
##  $ date  : chr  "2016_1" "2016_2" "2016_3" "2016_4" ...
##  $ aircon: int  434 504 769 420 759 1470 1542 651 469 336 ...
##  $ ice   : int  464 397 493 617 890 883 1292 1387 843 621 ...
##  $ kion  : num  5 4.9 8.3 14.1 19.4 22.4 26.2 27.5 23.8 18 ...
##  $ situdo: int  74 71 63 64 64 71 74 71 78 68 ...
head(x)
date aircon ice kion situdo
2016_1 434 464 5.0 74
2016_2 504 397 4.9 71
2016_3 769 493 8.3 63
2016_4 420 617 14.1 64
2016_5 759 890 19.4 64
2016_6 1470 883 22.4 71

売上の変化を図示します。

par(mfrow = c(2, 1))
plot(x$aircon, type='l',ylim=c(0,1700))
plot(x$ice, type='l' ,ylim=c(0,1700))

問１

airconとiceのグラフを見ると形は似通っていますが、iceのピークが少し右にずれているようです。
airconとiceは気温の影響を受けます。この場合、airconとiceには気温による擬似相関が発生しますが、このことは問５に出てきます。
iceの売上は、その日の気温に即時的に反応し、airconは前もって準備するので
airconのピークが少し前に来ると考えられます。

時系列分解

分析のため、データを時系列データに変換します。

aircon<-ts(x[,2],start=c(2016,1),frequency=12)
ice<-ts(x[,3],start=c(2016,1),frequency=12)
aircon
##       Jan  Feb  Mar  Apr  May  Jun  Jul  Aug  Sep  Oct  Nov  Dec
## 2016  434  504  769  420  759 1470 1542  651  469  336  427  571
## 2017  520  547  816  482  882 1560 1591  714  495  317  435  565
## 2018  515  582  916  534  901 1568 1763  818  557  376  493  627
## 2019  570  608 1000  651 1078 1662 1227  960  731  324  419  585
## 2020  481  551  905  598  947 1844 1464  940  632  379  494  635
ice
##       Jan  Feb  Mar  Apr  May  Jun  Jul  Aug  Sep  Oct  Nov  Dec
## 2016  464  397  493  617  890  883 1292 1387  843  621  459  562
## 2017  489  382  472  624  915  914 1394 1370  826  599  489  573
## 2018  507  416  607  746  894 1021 1506 1443  861  640  492  537
## 2019  494  423  542  667 1000  991 1236 1513  996  724  531  584
## 2020  510  482  610  689 1040 1123 1155 1658 1025  649  573  599

時系列データを図示すると、横軸は Time になります。

par(mfrow = c(2, 1))
ts.plot(aircon)
ts.plot(ice)

一般に、時系列データは、次の成分に分解されます。

観測値 (Obserbed) 実際の値
トレンド (Trend) データの長期的変化・傾向
季節変動 (Seasonal) 周期変動
誤差変動 (Random) 誤差による変化

トレンドは１２ヶ月ごとの移動平均で算出されます。
関数を利用して時系列分解（time series decomposition）を行います。

ts_decompose(ice)
ts_decompose(aircon)

トレンドを見ると、
アイスクリームは順調に売上を伸ばしているようですが、エアコンは2020にやや落ち込んでいるように見えます。
誤差変動は2019末に大きくマイナスに振れており、何か異常事態が起きていることが察せられます。ちょうどコロナの始まりと重なっているようです。

移動平均を具体的に見てみましょう。
トレンドは12ヶ月移動平均なので、SMA 関数を用いて12ヶ月移動平均に変換します。

i.ice<-SMA(ice,12)
i.ice
##           Jan      Feb      Mar      Apr      May      Jun      Jul      Aug
## 2016       NA       NA       NA       NA       NA       NA       NA       NA
## 2017 744.4167 743.1667 741.4167 742.0000 744.0833 746.6667 755.1667 753.7500
## 2018 755.4167 758.2500 769.5000 779.6667 777.9167 786.8333 796.1667 802.2500
## 2019 804.7500 805.3333 799.9167 793.3333 802.1667 799.6667 777.1667 783.0000
## 2020 809.7500 814.6667 820.3333 822.1667 825.5000 836.5000 829.7500 841.8333
##           Sep      Oct      Nov      Dec
## 2016       NA       NA       NA 742.3333
## 2017 752.3333 750.5000 753.0000 753.9167
## 2018 805.1667 808.5833 808.8333 805.8333
## 2019 794.2500 801.2500 804.5000 808.4167
## 2020 844.2500 838.0000 841.5000 842.7500
i.aircon<-SMA(aircon,12)
i.aircon
##           Jan      Feb      Mar      Apr      May      Jun      Jul      Aug
## 2016       NA       NA       NA       NA       NA       NA       NA       NA
## 2017 703.1667 706.7500 710.6667 715.8333 726.0833 733.5833 737.6667 742.9167
## 2018 743.2500 746.1667 754.5000 758.8333 760.4167 761.0833 775.4167 784.0833
## 2019 808.7500 810.9167 817.9167 827.6667 842.4167 850.2500 805.5833 817.4167
## 2020 810.5000 805.7500 797.8333 793.4167 782.5000 797.6667 817.4167 815.7500
##           Sep      Oct      Nov      Dec
## 2016       NA       NA       NA 696.0000
## 2017 745.0833 743.5000 744.1667 743.6667
## 2018 789.2500 794.1667 799.0000 804.1667
## 2019 831.9167 827.5833 821.4167 817.9167
## 2020 807.5000 812.0833 818.3333 822.5000

移動平均なので、2016の12月になって初めて値が現れます。

par(mfrow = c(2, 1))
ts.plot(i.aircon)
ts.plot(i.ice)

12ヶ月移動平均のグラフはトレンドと一致しています。

問２

移動平均がトレンドであることを知らないと、イメージしにくい問題です。
移動期間をいろいろ変えて図示してみます。
aircon

aircon9 <-SMA(aircon,9)
aircon6 <-SMA(aircon,6)
aircon3 <-SMA(aircon,3)
par(mfrow = c(3, 1))
ts.plot(aircon9)
ts.plot(aircon6)
ts.plot(aircon3)

⓪は３ヶ月移動平均　、②は９ヶ月移動平均　、③は６ヶ月移動平均
①は12ヶ月移動平均でトレンドになります。

問３

自己相関係数
未来は過去の影響を受けます。 時系列のある時刻のデータとkだけ離れた時刻のデータとの相関係数を、ラグkの自己相関係数と言い、\(\rho_k\) と書きます。
（数理的な説明は別の機会に譲ります）

自己相関(ACF)と偏自己相関(PACF)

ts_cor(aircon, lag.max = 60)

自己相関関数(Autocorrelation Function)は、１年ごとにピーク（赤）を迎え、時間の経過とともに小さくなっていきます。
偏自己相関(Partial AutoCorrelation Function)は、ラグhに関して時点t-hと時点tの間に存在するh-1個の観測値の影響を除去したあとの相関ですが、難しいのでここでは立ち入りません。 問３のグラフは自己相関(ACF)のグラフです。

ラグごとの散布図

ts_lags(aircon, lags = 1:12)

ラグ１で正の相関であったものが、ラグ５では負の相関になり、ラグ12（１年後）では強い正の相関となっており、これがACFのグラフのはじめの赤い棒に相当します。
この変化は、12ヶ月ごとに繰り返されます。

問４

相互相関(cross-correlation function)CCF
過去の自分との類似性を見るのが「自己相関」、他の時系列データとの類似性を見るのが「相互相関」です。
セールスアナリティクス様「時系列データの自己相関と相互相関をPythonで求めてみよう」参照

ccf(ice, aircon)