full_join()関数を利用して，株式データannual_stock_dataと財務データfinancial_dataを結合し，年次のデータフレームannual_dataを作成してみよう．
pacman::p_load(tidyverse)
# 財務データの読み込み
financial_data <- read_csv("../simulation_data/ch04_output.csv")
# 先週保存したannual_stock_data.csvを読み込み
annual_stock_data <- read_csv("../simulation_data/annual_stock_data.csv")
# 両データフレームの結合
annual_data <- annual_stock_data %>%
full_join(financial_data, by = c("firm_ID", "year")) # firm_IDとyearのペアをキーとして設定# 月次データの読み込み
stock_data <- read_csv("../simulation_data/ch05_output1.csv")
annual_data <- stock_data %>%
filter(month == 12) %>% # 12月のデータのみを抽出
select(year, firm_ID, ME) %>% # 追加したい列のみ選択
full_join(annual_data, ., by = c("year", "firm_ID")) %>% # 年次データと結合
mutate(ME = ME / 1e6) # 時価総額の単位を百万円に.) は，パイプ演算子%>%を用いてデータフレームを第一引数以外に代入する場合に用いられる．.) を用いてパイプ演算子から受け取ったデータフレームを第二引数に代入して，annual_dataをfull_join()関数 の第一引数として扱っている．geom_point()関数を用いるが，aes()関数の中にsize引数を追加すれば良いだけである．その関数内のalphaは透過度を表す引数であり，0から1までの適当な値を指定しよう．ここでは，仮に0.4を指定している．scale_size()関数を用いて時価総額を表す円の幅を指定している．annual_data %>%
filter(year == 2015,
firm_ID %in% 2:20, # firm_IDが2から20のデータを抽出
X > 0) %>% # 対数を取るため当期純利益（X）が正のデータのみ抽出
ggplot() +
geom_point(aes(x = log(sales), y = log(X), size = ME), alpha = 0.4) + # バブルチャートを描くにはsize引数を指定
scale_size(range = c(1, 20), name = "Market Equity") + # range引数でバブルの最小・最大面積を指定
scale_x_continuous(limits = c(8, 14)) + # 両軸の範囲を指定
scale_y_continuous(limits = c(2, 11)) +
theme_classic()
以上の仮定を満たす金融市場のことを，一般に完全資本市場（完全市場）と呼ぶ．「取引を行う上で完全に摩擦のない市場」といったイメージである．
\[ \begin{align*} \underbrace{\mathbb{E}\left[R_i\right]-R_F}_{\textbf{証券$i$のリスクプレミアム}} & =\underbrace{\beta_i}_{\textbf{証券$i$のマーケット・ベータ}}\underbrace{\left(\mathbb{E}\left[R_M\right]-R_F\right)}_{\textbf{市場リスクプレミアム}}\\ \text{ただし，}\quad \beta_i & =\frac{{\rm Cov}[R_i,\ R_M]}{{\rm Var}[R_M]} \end{align*} \]
\[ \begin{align} \underbrace{R_{i,t}}_{\textbf{証券$i$の実現リターン}} - \underbrace{R_{F,t}}_{\textbf{無リスク金利}} & = \beta_i(\underbrace{R_{M,t}}_{\substack{\textbf{市場ポートフォリオの} \\ \textbf{実現リターン}}} - \underbrace{R_{F,t}}_{\textbf{無リスク金利}}) +\varepsilon_{i,t}\nonumber\\ \underbrace{R_{i,t}^e}_{\textbf{証券$i$の実現超過リターン}} & = \beta_i\underbrace{R_{M,t}^e}_{\substack{\textbf{市場ポートフォリオの} \\ \textbf{実現超過リターン}}} +\varepsilon_{i,t} \label{eq:CAPM1} \tag{6.1} \end{align} \]
ここで，\(\varepsilon_{i,t}\)に関して以下の仮定を置こう．
lagged_BEMEを定義し，上記の基準に沿うデータだけを抽出し，firm_ID, year, Re, lagged_BEMEの四変数のみを残し，そのデータフレームをlm_sample_dataとして定義しよう．lagged_BEMEと実現超過リターンReの関係を散布図により可視化してみよう．lm_sample_data <- annual_data %>%
group_by(firm_ID) %>%
mutate(lagged_BEME = lagged_BE / lag(ME)) %>% # 前年度の簿価時価比率を計算
ungroup() %>%
filter(year == 2016, # firm_IDが1から10までの企業の2016年のデータを抽出
firm_ID <= 10) %>%
select(firm_ID, year, Re, lagged_BEME) %>% # 必要な列のみ抽出
drop_na() # 欠損データを削除
ggplot(lm_sample_data) +
geom_point(aes(x = lagged_BEME, y = Re)) + # x軸に簿価時価比率, y軸に超過リターン
labs(x = "BE/ME at the End of Year t", y = "Excess Return for Year t + 1") +
theme_classic()
geom_smooth()関数を用いて，散布図上に新しく直線を追加している．この関数の引数を見てみると，method引数で指定されている"lm"はlinear model（線形モデル）の略で，以下で詳細に説明する線形回帰を用いて直線を描くよう指定している．ggplot(lm_sample_data) +
geom_point(aes(x = lagged_BEME, y = Re)) +
geom_smooth(aes(x = lagged_BEME, y = Re),
method = "lm", se = FALSE, color = "black") + # 回帰直線を追加するにはgeom_smooth()関数を用いる
labs(x = "BE/ME at the End of Year t", y = "Excess Return for Year t + 1") +
theme_classic()
翌年のリターンと年度末のBE/MEを\(y\)と\(x\)に見立て，両者の関係を直線で表そうとした．この直線のことを線形モデルと呼ぶ．
そのモデルこそが，上式であり，モデルによって説明される変数\(y\)のことを被説明変数（従属変数）と呼び，他方，説明に利用する変数\(x\)のことを説明変数（独立変数）と呼ぶ．
被説明変数の平均的な動きの説明を試みようとする線形回帰モデルでは，残差\(e_{i}\)を不規則な動きをする確率変数の実現値として捉える．そして，その確率変数のことを誤差項 (error term)とか，単に誤差と呼び，残差と区別するために\(\varepsilon_{i}\)として定義しよう．一般に線形回帰モデルでは，誤差項をうまくモデルに取り入れて，\((X_i, Y_i)\)の間に以下のような関係が成立すると想定する． \[ \begin{align} Y_i = \beta_0 + \beta_1 X_i + \varepsilon_i \label{eq:linear_model} \tag{5.5} \end{align} \] 誤差項\(\varepsilon_{i}\)の条件付き期待値\(\mathbb{E}[\varepsilon_{i} \mid X_{i}]\)はゼロと仮定する．
ここで取り上げた例だと，\(X_i\)は\(t\)年度末の簿価時価比率\(\mathit{BE/ME}_{i,t}\)で，\(Y_i\)は\(t+1\)年の無リスク金利に対する超過リターン\(R_{i,t+1}^e=R_{i,t+1} - R_{F,t+1}\)，添字の\(i\)はfirm_IDに該当する．したがって，より具体的に書けば， \[
\begin{align*}
\underbrace{R_{i,t+1}^e}_{Y_{i}} = \beta_{0} + \beta_{1} \underbrace{\mathit{BE/ME}_{i,t}}_{X_{i}} + \underbrace{\varepsilon_{i,t+1}}_{\varepsilon_{i}}
\end{align*}
\]
(\(\ref{eq:linear_model}\))式において，\((\beta_0, \beta_1)\)は回帰係数と呼ばれ，\(X_i\)と\(Y_i\)の線形関係を表すパラメータである．
回帰係数の値はデータから直接観察することはできないので，線形回帰モデルではその推定値が関心の対象となる．
例えば，先ほどの散布図では，\(\hat{\beta}_0 \approx 0.055\)，\(\hat{\beta}_1 \approx 0.076\)であるが，ここで\((\hat{\beta}_0, \hat{\beta}_1)\)という変数名を使ったのは，\((\beta_0, \beta_1)\)の推定値であることを強調するためである．
観察されたデータから\((\beta_0, \beta_1)\)を推定する方法として最もよく用いられるのが最小二乗法 (Ordinary Least Squares)，略してOLSである．
OLSでは，誤差の二乗和を最小にするように，観察されたデータから回帰係数を推定する．観測値\(i\)の誤差は(\(\ref{eq:linear_model}\))式より\(\varepsilon_i = Y_i - (\beta_0 + \beta_1 X_i)\)として求めることができるから，データ数が\(n\)ならば，誤差の二乗和は， \[\begin{align*} \sum_{i=1}^n \varepsilon_i^{2} = \sum_{i=1}^n \left\{ Y_i - ({\beta_0} + {\beta_1} X_i ) \right\}^2 \end{align*}\]
となり，\(\beta_{0}\)と\(\beta_{1}\)の関数で表せることが分かる．OLSでは， \[\begin{align*} \min_{\hat{\beta}_0, \hat{\beta}_1} \sum_{i=1}^n \left\{ Y_i - (\hat{\beta}_0 + \hat{\beta}_1 X_i ) \right\}^2 \end{align*}\]
を解く最小二乗推定量\((\hat{\beta}_0, \hat{\beta}_{1})\)によって，未知の回帰係数\((\beta_0, \beta_{1})\)を推定する．
lm()関数では，第一引数で~演算子を用いて回帰モデルを記述する．data引数でそのデータフレームそのものを指定すると，回帰モデルではデータフレーム名は省略して，列名のみの記述が可能になる． lm()関数は，返り値として名前付きのリストを返すので，以下ではnames()関数を用いて，その一覧を表示させている．lm_results <- lm(Re ~ lagged_BEME, data = lm_sample_data) # ~の左に従属変数, 右に独立変数を記す
names(lm_results) [1] "coefficients" "residuals" "effects" "rank"
[5] "fitted.values" "assign" "qr" "df.residual"
[9] "xlevels" "call" "terms" "model" lm()関数で計算した回帰係数を確認するためには，リストlm_results内の一つ目の要素をprint()関数により表示すれば良い．リストの要素へのアクセスの仕方は，第2回資料を参照．print(lm_results[[1]])(Intercept) lagged_BEME
0.05513385 0.07552921 $演算子を使って以下のようにアクセスすることも可能．print(lm_results$coefficients)(Intercept) lagged_BEME
0.05513385 0.07552921 broomパッケージを利用した結果の確認方法lm()関数は線形回帰の結果をリスト形式で詳細に返す．ただし，返り値がリストのままだと分析が不便であるので，代わりに外部パッケージのbroomを用いて，データフレームに変換しよう．tidy()関数が準備されている．# broomパッケージの読み込み
pacman::p_load(broom)
# lm_resultsをデータフレーム形式に変換
tidy(lm_results) | モデル | \(\underset{(t\text{値})}{\hat{\beta}_{1}}\) | 解釈 |
|---|---|---|
| \(Y_{i} = \beta_{0} + \beta_{1}X_{i} + \varepsilon_{i}\) | \(\underset{(4.33)}{32.6}\) | 今期のR&Dが1百万円上昇すれば，翌期の売上高は平均的に32.6百万円上昇する |
| \(Y_{i} = \beta_{0} + \beta_{1}\log(X_{i}) + \varepsilon_{i}\) | \(\underset{(6.40)}{36,129}\) | 今期のR&Dが1%上昇すれば，翌期の売上高は平均的に361.29百万円上昇する |
| \(\log(Y_{i}) = \beta_{0} + \beta_{1}X_{i} + \varepsilon_{i}\) | \(\underset{(10.50)}{0.00055}\) | 今期のR&Dが1百万円上昇すれば，翌期の売上高は平均的に0.055%上昇する |
| \(\log(Y_{i}) = \beta_{0} + \beta_{1}\log(X_{i}) + \varepsilon_{i}\) | \(\underset{(11.62)}{0.497}\) | 今期のR&Dが1%上昇すれば，翌期の売上高は平均的に0.497%上昇する |
株主資本や利益など，企業のファンダメンタルズを表す数値をベンチマークとして株価（延いては時価総額）が割安な銘柄をバリュー株と呼び，その対極にあるのはグロース株と呼ぶ．↩︎