人間自由度指数は、世界中の多くの国々の様々な変数を通じて「自由」という考え方を要約しようとするレポートです。 それは、政治的、宗教的、経済的、個人的な自由といった異なる種類の自由と、その他の社会的、経済的状況との関係を大まかに客観的に測定する役割を果たします。 人間自由度指数は、カトー研究所、フレーザー研究所、フリードリヒ・ナウマン財団の自由主義研究所によって毎年共同出版されるレポートです。

このラボでは、人間自由度指数レポートのデータを分析します。 あなたの目的は、データ内のいくつかの関係をグラフィカルおよび数値的に要約し、自由に関する物語を語るのに役立つ変数を特定することです。

はじめに

パッケージの読み込み

このラボでは、tidyverse スイートのパッケージを使用してデータを探索および視覚化します。 また、残差平方和を最小化する回帰直線を特定するためにstatsrパッケージを使用し、回帰出力を整理するためにbroomパッケージを使用します。 データはOpenIntroリソースのコンパニオンパッケージであるopenintroパッケージで見つけることができます。

パッケージを読み込みましょう。

library(tidyverse)
library(openintro)
library(statsr)
library(broom)

再現可能なラボレポートの作成

新しいラボレポートを作成するには、RStudioで「New File」->「R Markdown…」に移動し、「From Template」を選択してから、テンプレートのリストから「Lab Report for OpenIntro Statistics Labs」を選択します。

データ

私たちが使用するデータはopenintroパッケージにあり、人間自由度指数を表すhfiと呼ばれます。

  1. データセットの次元は何ですか？ 各行は何を表していますか？

  2. データセットは多くの年にまたがっていますが、2016年のデータにのみ関心があります。 hfiデータフレームを2016年にフィルタリングし、6つの変数を選択して、結果をhfi_2016という名前のデータフレームに割り当ててください。

  3. 個人の自由度スコアpf_scoreとpf_expression_controlの関係を表示するために、どのような種類のプロットを使用しますか？pf_expression_controlを変数として使用してこの関係をプロットしてください。関係は線形に見えますか？もしある国のpf_expression_control（メディア内容への政治的圧力と統制のスコアで、10点満点中0が最も高い）を知っていたら、個人自由度スコアを予測するために線形モデルを使用しても問題ないと思いますか？

関係が線形に見える場合、相関係数を使用して関係の強さを定量化できます。

hfi_2016 |>
  summarise(cor(pf_expression_control, pf_score))

残差平方和

このセクションでは、「残差平方和」とは何を意味するかを調査するためにインタラクティブな関数を使用します。 この関数は、RStudioのコンソールで実行する必要があり、R Markdownドキュメントでは実行できません。 この関数を実行するには、hfiデータセットが環境にロードされている必要もあります。 また、右下隅の「Plots」タブにプロットを作成するのに十分なスペースがあることを確認する必要があります。

単一変数の分布を記述した方法を思い出してください。 中心、広がり、形状といった特性について議論しました。 pf_expression_controlとpf_scoreのような2つの数値変数の関係を記述できることも有用です。

  1. 前の演習のプロットを見て、これら2つの変数の関係を記述してください。関係の形式、方向、強さ、および異常な観測について議論するようにしてください。

単一変数を要約するために平均と標準偏差を使用したのと同様に、関連性を最もよく表す直線を見つけることで、これら2つの変数の関係を要約できます。 以下の対話型関数を使用して、点の雲を最もよく通過すると考える線を選択してください。

plot_ss(x = pf_expression_control, y = pf_score, data = hfi_2016)

このコマンドを実行すると、プロット上に2つの点をクリックして線定義するよう求められます。 それが完了すると、指定した線が黒で、残差が青で表示されます。

プロットがコードチャンクの下に表示され、線を作成するために点を選択できない場合は、次の手順を実行してください。

  • RStudioの上部にあるツールバーに移動します。
  • 「Global Options…」をクリックします。
  • 「R Markdown pane」（左側）をクリックします。
  • 「Show output inline for all R Markdown documents」というボックスのチェックを外します。

残差は、観測値と線によって予測された値との差であることを思い出してください。

\[ e_i = y_i - \hat{y}_i \]

線形回帰を行う最も一般的な方法は、残差平方和を最小化する線を選択することです。 平方残差を視覚化するには、プロットコマンドを再実行し、引数showSquares = TRUEを追加します。

plot_ss(x = pf_expression_control, y = pf_score, data = hfi_2016, showSquares = TRUE)

plot_ss関数の出力は、線の傾きと切片、および平方和を提供することに注意してください。

  1. plot_ssを使用して、平方和を最小限に抑えるのに適した線を選択してください。関数を数回実行してください。得られた最小の平方和はいくつでしたか？隣のチームと比較してどうでしたか？

線形モデル

試行錯誤によって正しい最小二乗線、すなわち残差平方和を最小化する線を得ようとすることは、かなり面倒です。 代わりに、Rのlm関数を使用して線形モデル（別名回帰直線）を適合させることができます。

m1 <- lm(pf_score ~ pf_expression_control, data = hfi_2016)

lm()関数の最初の引数は、y ~ xの形式を取る式です。 ここでは、pf_scoreをpf_expression_controlの関数として線形モデルを作成したいと読むことができます。 2番目の引数は、Rが2つの変数を見つけるためにhfiデータフレーム内を検索することを指定します。

注: データフレームが最初の引数ではないため、パイピングはここでは機能しません！

lm()の出力は、適合したばかりの線形モデルに関するすべての情報を含むオブジェクトです。 tidy()関数を使用してこの情報にアクセスできます。

tidy(m1)

この出力を部分ごとに見ていきましょう。 まず、モデルを記述するために使用された式が、「Call」として表示される上部に示されています。 式の後には、残差の5数要約があります。 次に示される「Coefficients」テーブルが重要です。その最初の列には、線形モデルのy切片とpf_expression_controlの係数が表示されます。 このテーブルを使用して、線形モデルの最小二乗回帰直線を記述できます。

\[ \hat{y} = 4.28 + 0.542 \times pf\_expression\_control \]

この方程式は2つのことを教えてくれます。

  • pf_expression_controlが0の国（メディア内容への政治的圧力が最も大きい国）の場合、平均個人自由度スコアは4.28になると予想されます。
  • pf_expression_controlが1単位増加するごとに、国の平均個人自由度スコアが0.542単位増加すると予想されます。

説明変数によって説明される応答変数の変動の割合である\(R^2\)を使用してモデルの適合度を評価できます。 この情報にアクセスするにはglance()関数を使用します。

glance(m1)

このモデルの場合、pf_scoreの変動の71.4%はpf_expression_controlによって説明されます。

  1. pf_expression_controlを使用してhf_score、つまり合計人間自由度スコアを予測する新しいモデルを適合させてください。Rの出力からの推定値を使用して、回帰直線の方程式を記述してください。人間自由度とメディアコンテンツへの政治的圧力の間の関係の文脈で、傾きは何を教えてくれますか？

予測と予測誤差

m1の最小二乗線を重ねた散布図を作成しましょう。

ggplot(data = hfi_2016, aes(x = pf_expression_control, y = pf_score)) +
  geom_point() +
  geom_smooth(method = "lm", se = FALSE)

ここでは、文字通りプロットの上にレイヤーを追加しています。 geom_smoothは線形モデルを適合させることで線を作成します。 また、線に関連する標準誤差seを表示することもできますが、ここではそれを抑制します。

この線は、任意の\(x\)値における\(y\)を予測するために使用できます。 観測データの範囲を超える\(x\)値に対して予測を行うことを外挿と呼び、通常は推奨されません。 しかし、データ範囲内でなされる予測はより信頼性があります。 それらは残差を計算するためにも使用されます。

  1. もし誰かが実際のデータではなく最小二乗回帰線を見た場合、pf_expression_controlの評価が3の国の個人自由度スコアをどのように予測するでしょうか？これは過大評価ですか、過小評価ですか、またどの程度ですか？言い換えれば、この予測の残差は何ですか？

モデル診断

線形モデルが信頼できるかどうかを評価するには、(1)線形性、(2)ほぼ正規の残差、(3)一定の変動性の3つを確認する必要があります。

これらのチェックを行うためには、適合値（予測値）と残差にアクセスできる必要があります。 これらを計算するにはaugment()関数を使用できます。

m1_aug <- augment(m1)

線形性: pf_scoreとpf_expression_controlの関係が線形であるかどうかは、散布図を使用してすでに確認しました。 残差対適合値（予測値）のプロットでもこの条件を確認する必要があります。

ggplot(data = m1_aug, aes(x = .fitted, y = .resid)) +
  geom_point() +
  geom_hline(yintercept = 0, linetype = "dashed", color = "red") +
  xlab("Fitted values") +
  ylab("Residuals")

ここで注意すべきは、m1もデータセットとして機能できることです。なぜなら、その中に適合値（\(\\\hat{y}\)）と残差が保存されているからです。 また、ここではコードを凝っていることにも注意してください。 最初のレイヤーで散布図を作成した後（コードの最初の行）、y = 0に赤い水平破線（残差が0の周りに分布しているかを確認するのに役立つ）を重ね、軸ラベルもより分かりやすくするために名前を変更しています。

  1. 残差プロットに何らかの明らかなパターンはありますか？これは2つの変数の関係の線形性について何を示していますか？


ほぼ正規の残差: この条件を確認するには、残差のヒストグラムを見ることができます。

ggplot(data = m1_aug, aes(x = .resid)) +
  geom_histogram(binwidth = 0.25) +
  xlab("Residuals")
  1. ヒストグラムに基づいて、ほぼ正規の残差という条件は破られているように見えますか？なぜそう思いますか？


一定の変動性:

  1. 残差対適合プロットに基づいて、一定の変動性という条件は破られているように見えますか？なぜそう思いますか？

さらなる演習

  • pf_scoreと強く相関すると考える別の変数を選択してください。 2つの変数の散布図を作成し、線形モデルを適合させてください。 一見して、線形関係があるように見えますか？

  • この関係はpf_scoreとpf_expression_controlの関係とどのように比較されますか？ 2つのモデル要約の\(R^2\)値を使用して比較してください。 あなたの独立変数はpf_scoreをよりよく予測しているように見えますか？ なぜそう思いますか？

  • 適切な視覚化を使用してモデル診断を確認し、モデル条件が満たされているかどうかを評価してください。

  • 関心のある別の変数ペアを選択し、それらの関係を視覚化してください。 その関係に驚きましたか、それとも予想通りでしたか？ なぜこれらの変数に関心を持ったのか、そして観察された関係に驚いた（または驚かなかった）理由について議論してください。

Creative Commons License
この作品は、クリエイティブ・コモンズ 表示 - 継承 4.0 国際ライセンスの下に提供されています。

```