7 第7回 回帰分析で新店舗の売上予測
7.1 準備
ここでは週次集計POSデータを使って、他の店で売れている商品を見つけ出す方法を学びます。 まずパッケージとデータを準備します。
第7回ファイルで使うデータはchp7.xlsxです。 いつものようにchp7.xlsxのシート名を確認します。
readxl::excel_sheets("data/chp7.xlsx") [1] "2023年11月期データ" "図7-1-図7-5" "図7-7"
[4] "図7-11-図7-15" "図7-18-図7-20" "図7-22-図7-25"
[7] "図7-27" "図7-28" "表7-2"
[10] "図7-29" "図7-30" "表7-3"
[13] "図7-31" "表7-4" "図7-32"
[16] "表7-5" "表7-6" "図7-33"
[19] "図7-35" たくさんシートがあるMS Excelファイルですが、1番目の「2023年11月期データ」のシートを読み込みます。 readxlパッケージのread_excel()関数にsheet = 1を指定することで、1番目のシートを読み込むことができます。
df <- readxl::read_excel("data/chp7.xlsx", sheet = 1)
head(df) # 先頭6行を表示| 店舗 | 店舗面積 | 商圏人口 | 最寄り駅からの距離 | 駐車場台数 | 競合店舗数 | 月間売上高 |
|---|---|---|---|---|---|---|
| A001 | 1500 | 207.6 | 100 | 27 | 12 | 130069 |
| A002 | 1010 | 210.8 | 590 | 20 | 3 | 72981 |
| A003 | 1260 | 174.5 | 180 | 57 | 19 | 84050 |
| A004 | 1180 | 288.8 | 580 | 57 | 9 | 108947 |
| A005 | 1160 | 216.1 | 1170 | 65 | 16 | 65733 |
| A006 | 1290 | 251.0 | 570 | 43 | 2 | 114216 |
このデータフレームには、
- 店舗 : 文字列
- 店舗面積 : 数値
- 商圏人口 : 数値
- 最寄り駅からの距離 : 数値
- 駐車場台数 : 数値
- 競合店舗数 : 数値
- 月間売上高 : 数値
という7つの変数が含まれています。
それぞれの詳細はテキストp.215を参照してください。
7.2 散布図
2変数の関係を視覚的に表すグラフとして散布図(scatter diagram)があります。 サンプルサイズが \(i = 1,\dots, n\) この変数 \(X\) と \(Y\) の\(i\)番目の観測値を \((x_i, y_i)\) とすると、散布図はこれらの点を座標平面上にプロットしたものです。
まず最寄り駅からの距離と月間売上高の散布図を描いてみます。 フォントの設定や見た目の設定をmystyleとしてlistにまとめておきます。
mystyle <- list(
theme_calc(base_family = "HiraKakuPro-W3"), # フォントの設定
scale_fill_tableau(name = "Tableau 20") # カラーパレットの設定
)ggplot2パッケージのggplot()関数を使い、geom_point()で散布図を描きます。
df |>
ggplot() + aes(x = 最寄り駅からの距離, y = 月間売上高) + # 軸の設定
geom_point() + # 散布図の描画
labs(title = "図7-1 最寄り駅からの距離と月間売上高の散布図") + mystyle
いろんな変数の組み合わせの散布図を描いてみます。 ここでは教科書に合わせて原点を0にしていますが、散布図ではかならずしも原点を0にする必要なないので、軸の範囲を指定しなくてもよいです。
g1 <- df |>
ggplot() + aes(x = 店舗面積, y = 月間売上高) +
geom_point() +
xlim(0, 2000) + ylim(0,160000) +
labs(title = "図7-2 店舗面積と月間売上高") + mystyle
g2 <- df |>
ggplot() + aes(x = 商圏人口, y = 月間売上高) +
geom_point() +
xlim(0, 350) + ylim(0,160000) +
labs(title = "図7-3 商圏人口と月間売上高") + mystyle
g3 <- df |>
ggplot() + aes(x = 駐車場台数, y = 月間売上高) +
geom_point() +
xlim(0, 90) + ylim(0,160000) +
labs(title = "図7-4 駐車場台数と月間売上高") + mystyle
g4 <- df |>
ggplot() + aes(x = 競合店舗数, y = 月間売上高) +
geom_point() +
xlim(0, 25) + ylim(0,160000) +
labs(title = "図7-5 競合店舗数と月間売上高") + mystyle
(g1 + g2) / (g3 + g4)
7.3 相関係数
2変数間の関係を数値で表す指標として相関係数(correlation coefficient)があります。
相関係数という尺度にはいろんな種類があり、最もよく利用されているのがピアソンの積率相関係数です。 他にも、スピアマンの順位相関係数がありますが、最もよく利用されるのがピアソンの積率相関係数なので、 ここでは相関係数といえばピアソンの積率相関係数を指すことにします。
店舗面積、商圏人口、駐車場台数、競合店舗数と月間売上高の相関係数を計算してみます。 Rには相関係数を返す基本関数としてcor()があります。 cor()関数とcorrplotパッケージのcorrplot()関数を使って相関係数行列を可視化します。
par(family= "HiraKakuProN-W3") # macの文字化け対策
df |>
select(店舗面積, 商圏人口, 最寄り駅からの距離, 駐車場台数, 競合店舗数, 月間売上高) |> # 変数を選択
cor() |> # 相関係数行列を計算
corrplot::corrplot(method = "number") # 相関係数行列を可視化
月間売上高と相関関係があるのは、最寄り駅からの距離(\(-0.71\))と、商圏人口(\(0.45\))といえます。
7.4 回帰分析
次に、月間売上高を説明するための回帰モデルを構築します。 回帰分析では、結果を表す変数と原因を表す変数の関係、つまり因果関係(causal relationship)を調べます。
| 種類 | 結果を表す変数・出力 | 原因を表す変数・入力 |
|---|---|---|
| パターン1 | 従属変数 (dependent variable) | 独立変数 (independent variables) |
| パターン2 | 被説明変数 (explaned variable) | 説明変数 (explanatory variables) |
| パターン3 | 目的変数 (objective variable) | 説明変数 (explanatory variables) |
| パターン4 | 応答変数 (response variable) | 予測変数 (predictor variable) |
用語法は文脈によって異なります。ここでの分類は完全に松浦の主観です。
7.4.1 回帰モデル
回帰モデルは、従属変数と独立変数の関係を数学的に表現するモデルを構築することです。 従属変数を \(1 \times n\) ベクトル\(y\)で、独立変数を \(k \times n\) 行列 \(X\)とすると、回帰モデルは以下のように表されます。
\[ \boldsymbol{y} = \boldsymbol{X} \boldsymbol{\beta} + \boldsymbol{\varepsilon} \]
ここで、\(\boldsymbol{\beta}\)は \(k \times 1\) ベクトルで、\(\boldsymbol{\varepsilon}\)は \(1 \times n\) ベクトルで、誤差項を表します。 成分を書き下すと以下のようになります。
\[ \begin{aligned} \begin{pmatrix} y_1 \\ y_2 \\ \vdots \\ y_n \end{pmatrix} = \begin{pmatrix} 1 & x_{11} & x_{12} & \cdots & x_{1k} \\ 1 & x_{21} & x_{22} & \cdots & x_{2k} \\ \vdots & \vdots & \vdots & \ddots & \vdots \\ 1 & x_{n1} & x_{n2} & \cdots & x_{nk} \\ \end{pmatrix} \begin{pmatrix} \beta_0 \\ \beta_1 \\ \vdots \\ \beta_k \end{pmatrix} + \begin{pmatrix} \varepsilon_1 \\ \varepsilon_2 \\ \vdots \\ \varepsilon_n \end{pmatrix} \end{aligned} \]
これを展開すると、よく教科書に出てくるような次の式になります。
\[ \begin{aligned} y_i & = \beta_0 + \beta_1 x_{i1} + \beta_2 x_{i2} + \cdots + \beta_k x_{ik} + \varepsilon_i \\ & = \beta_0 + \sum_{j=1}^{k} \beta_j x_{ij} + \varepsilon_i \end{aligned} \]
とりわけ\(k=1\)の場合、つまり独立変数が1つだけの回帰モデルを単回帰モデルと呼ばれます。 テキストでは単回帰を繰り返していますが、通常は複数の独立変数を組み込んだ重回帰モデルを構築します。
\[ \begin{aligned} \text{月間売上高}_i & = \beta_0 + \beta_1 店舗面積_i + \beta_2 商圏人口_i + \\ & + \beta_3 最寄り駅からの距離_i + \beta_4 駐車場台数_i + \\ & + \beta_5 競合店舗数_i + \varepsilon_i \end{aligned} \]
このモデルをRで推定します。
# 回帰モデルを`reg_model`に格納
reg_model <- "月間売上高 ~ 店舗面積 + 商圏人口 + 最寄り駅からの距離 + 駐車場台数 + 競合店舗数"
OLS <- lm(reg_model, data = df) # 最小二乗法で回帰分析
summary(OLS) # 結果を出力
Call:
lm(formula = reg_model, data = df)
Residuals:
Min 1Q Median 3Q Max
-22884.9 -8871.7 594.1 9150.0 24149.5
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) 70086.274 9100.157 7.702 0.000000000011489 ***
店舗面積 4.403 4.763 0.924 0.358
商圏人口 199.768 24.087 8.294 0.000000000000636 ***
最寄り駅からの距離 -47.798 3.691 -12.951 < 0.0000000000000002 ***
駐車場台数 9.441 54.552 0.173 0.863
競合店舗数 -159.811 195.141 -0.819 0.415
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Residual standard error: 11470 on 97 degrees of freedom
Multiple R-squared: 0.7144, Adjusted R-squared: 0.6997
F-statistic: 48.53 on 5 and 97 DF, p-value: < 0.00000000000000022Estimateは係数の推定値、Std. Errorは標準誤差、t valueはt値、Pr(>|t|)はp値です。 ここでは、商圏人口と最寄り駅からの距離が0.1％水準で統計的に有意となり、商圏人口が月間売上高に正の影響を、最寄り駅からの距離が負の影響を与えていることが分かります。 商圏人口は人数、最寄り駅からの距離は長さといったように変数ごとに単位が異なっているため、回帰係数の大小を比較することはできません。 しかしデータを標準化してから回帰分析を行うと、各変数の係数が標準偏差単位で解釈できるので、変数間の比較がしやすくなります。
7.4.2 標準化
標準化(standardization)は、変数の単位を揃えるために行います。 確率変数\(X\)を標準化するには、平均\(\bar{X}\)を引いて標準偏差\(s\)で割ります。
\[ \begin{aligned} Z_i & = \frac{X_i - \bar{X}}{s} \end{aligned} \]
この標準化された確率変数\(Z\)は平均が0、標準偏差が1になります。 したがう分布は元のままですが，平均が0，標準偏差が1になるため，異なる尺度の変数を比較する際に便利です。
月刊売上高で確認してみましょう。 月間売上高の平均は92501.56、標準偏差は20931.89です。 標準化を行うために便利なパッケージであるscale()関数を使います。
df <- df |>
mutate(
標準化月間売上高 = scale(月間売上高) # 月間売上高を標準化
)標準化月間売上高の平均は0、標準偏差は1となっていることが確認できました。 分布を図で表すためにヒストグラムを描いてみます。
hist(df$月間売上高, main = "月間売上高の分布", xlab = "月間売上高") # ヒストグラム
abline(v = mean(df$月間売上高), col = "red", lwd = 2) # 平均値に赤線を追加
正規化したヒストグラムを描いてみます。
hist(scale(df$月間売上高), main = "月間売上高の標準化", xlab = "月間売上高") # 標準化したヒストグラム
abline(v = mean(scale(df$月間売上高)), col = "red", lwd = 2) # 平均値に赤線を追加
正規化する前とした後の回帰分析の結果を比較してみましょう。 ここでは、modelsummaryパッケージを使って結果をまとめて表示します。
df_std <- df |>
select(店舗面積, 商圏人口, 最寄り駅からの距離, 駐車場台数, 競合店舗数, 月間売上高) |>
scale() |> # 標準化
as.data.frame() # データフレームに変換
result_std <- lm(月間売上高 ~ 店舗面積 + 商圏人口 + 最寄り駅からの距離 + 駐車場台数 + 競合店舗数, data = df_std)
results <- list("OLS" = OLS, "標準化" = result_std)
modelsummary::msummary(results,
stars = TRUE,
digits = 3)| OLS | 標準化 | |
|---|---|---|
| + p | ||
| (Intercept) | 70086.274*** | 0.000 |
| (9100.157) | (0.054) | |
| 店舗面積 | 4.403 | 0.052 |
| (4.763) | (0.056) | |
| 商圏人口 | 199.768*** | 0.465*** |
| (24.087) | (0.056) | |
| 最寄り駅からの距離 | -47.798*** | -0.715*** |
| (3.691) | (0.055) | |
| 駐車場台数 | 9.441 | 0.010 |
| (54.552) | (0.056) | |
| 競合店舗数 | -159.811 | -0.045 |
| (195.141) | (0.055) | |
| Num.Obs. | 103 | 103 |
| R2 | 0.714 | 0.714 |
| R2 Adj. | 0.700 | 0.700 |
| AIC | 2225.7 | 176.2 |
| BIC | 2244.2 | 194.7 |
| Log.Lik. | -1105.860 | -81.110 |
| F | 48.527 | 48.527 |
| RMSE | 11131.88 | 0.53 |