このラボでは、統計学で最も重要な確率分布である正規分布を調査します。 データがほぼ正規分布していると確信できれば、多くの強力な統計的手法への扉が開かれます。 ここでは、Rのグラフィカルツールを使用してデータの正規性を評価し、正規分布から乱数を生成する方法も学習します。

はじめに

パッケージの読み込み

このラボでは、tidyverse スイートのパッケージとopenintro パッケージを使用してデータを探索および視覚化します。

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

library(tidyverse)
library(openintro)

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

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

データ

今週はファストフードデータを扱います。 このデータセットには、世界で最も人気のあるファストフードレストランのメニュー項目515件のデータが含まれています。 データの最初の数行を簡単に見てみましょう。

以前のようにglimpseを使用するか、headを使用してこれを行うことができます。

library(tidyverse)
library(openintro)
head(fastfood)
## # A tibble: 6 × 17
##   restaurant item       calories cal_fat total_fat sat_fat trans_fat cholesterol
##   <chr>      <chr>         <dbl>   <dbl>     <dbl>   <dbl>     <dbl>       <dbl>
## 1 Mcdonalds  Artisan G…      380      60         7       2       0            95
## 2 Mcdonalds  Single Ba…      840     410        45      17       1.5         130
## 3 Mcdonalds  Double Ba…     1130     600        67      27       3           220
## 4 Mcdonalds  Grilled B…      750     280        31      10       0.5         155
## 5 Mcdonalds  Crispy Ba…      920     410        45      12       0.5         120
## 6 Mcdonalds  Big Mac         540     250        28      10       1            80
## # ℹ 9 more variables: sodium <dbl>, total_carb <dbl>, fiber <dbl>, sugar <dbl>,
## #   protein <dbl>, vit_a <dbl>, vit_c <dbl>, calcium <dbl>, salad <chr>

すべての観測値に対して17の測定値があり、その多くは栄養成分であることがわかります。

まず、レストラン、カロリー、脂肪からのカロリーの3つの列に焦点を当てます。

まず、マクドナルドとデイリークイーンの製品のみに焦点を当てましょう。

mcdonalds <- fastfood |>
  filter(restaurant == "Mcdonalds")
dairy_queen <- fastfood |>
  filter(restaurant == "Dairy Queen")
  1. これら2つのレストランのオプションにおける脂肪からのカロリー量の分布を視覚化するプロット(または複数のプロット)を作成してください。それらの中心、形状、広がりはどのように比較されますか?

正規分布

分布の説明で、「ベル型」や「正規」といった言葉を使用しましたか? 単峰性の対称分布に直面すると、そう言いたくなるものです。

その説明がどれほど正確であるかを確認するには、ヒストグラムの上に正規分布曲線をプロットして、データが正規分布にどれだけ密接に従っているかを確認できます。 この正規曲線は、データと同じ平均と標準偏差を持つ必要があります。 デイリークイーン製品の脂肪からのカロリーに焦点を当てるので、それらを別々のオブジェクトとして保存し、後で参照する統計量を計算しましょう。

dqmean <- mean(dairy_queen$cal_fat)
dqsd <- sd(dairy_queen$cal_fat)

次に、背景として使用する密度ヒストグラムを作成し、lines関数を使用して正規確率曲線を重ね合わせます。 度数ヒストグラムと密度ヒストグラムの違いは、度数ヒストグラムでは棒の高さが観測値の合計数に加算されるのに対し、密度ヒストグラムでは棒の面積が1に加算されることです。 各棒の面積は、単に棒の高さに幅を掛けたものとして計算できます。 密度ヒストグラムを使用すると、曲線も曲線下の面積が1の正規確率密度関数であるため、ヒストグラムの上に正規分布曲線を適切に重ね合わせることができます。 度数ヒストグラムと密度ヒストグラムはどちらもまったく同じ形状を示します。それらはy軸が異なるだけです。 これは、以前に作成した度数ヒストグラムと以下のコマンドで作成された密度ヒストグラムを比較することで確認できます。

ggplot(data = dairy_queen, aes(x = cal_fat)) +
  geom_blank() +
  geom_histogram(aes(y = ..density..)) +
  stat_function(fun = dnorm, args = c(mean = dqmean, sd = dqsd), col = "tomato")
## Warning: The dot-dot notation (`..density..`) was deprecated in ggplot2 3.4.0.
## ℹ Please use `after_stat(density)` instead.
## This warning is displayed once every 8 hours.
## Call `lifecycle::last_lifecycle_warnings()` to see where this warning was
## generated.

geom_blank()で空白のプロットを初期化した後、ggplot2パッケージ(tidyverse内)を使用すると、追加のレイヤーを追加できます。 最初のレイヤーは密度ヒストグラムです。 2番目のレイヤーは統計関数、つまり正規曲線の密度であるdnormです。 曲線が脂肪からのカロリーの列と同じ平均と標準偏差を持つように指定します。 col引数は、描画される線の色を単に設定します。 もし省略した場合、線は黒で描画されます。

  1. このプロットに基づいて、データはほぼ正規分布に従っているように見えますか?

正規分布の評価

ヒストグラムの形状を目で見てデータがほぼ正規分布しているかどうかを判断する方法は1つですが、ヒストグラムが曲線にどれだけ近いかを判断するのは困難な場合があります。 別の方法として、正規確率プロット、別名「分位点-分位点」の正規Q-Qプロットを作成する方法があります。

ggplot(data = dairy_queen, aes(sample = cal_fat)) +
  geom_line(stat = "qq")

今回は、geom_line()レイヤーを使用し、stat引数でQ-Qプロットを作成することを指定します。 ここで重要なのは、aes()内でxを使用する代わりにsampleを使用する必要があることです。

x軸の値は、平均0、標準偏差1の理論的な正規曲線(すなわち、標準正規分布)の分位点に対応します。 y軸の値は、元の標準化されていないサンプルデータの分位点に対応します。 しかし、サンプルデータを標準化したとしても、Q-Qプロットは同じに見えます。 ほぼ正規分布しているデータセットは、点が対角線に密接に沿う確率プロットになります。 正規性からの逸脱は、これらの点がその線から逸脱することにつながります。

デイリークイーンの脂肪からのカロリーのプロットは、線に沿う傾向があるが、上端に向かっていくつかの逸脱した点があることを示しています。 上記ヒストグラムで遭遇したのと同じ問題に直面します。どれだけ近ければ十分なのか?

この質問に対処する有用な方法は、それを「私が正規分布から得たと分かっているデータの確率プロットはどのように見えるか?」と再フレーズすることです。 これは、rnormを使用して正規分布からデータをシミュレートすることで答えることができます。

sim_norm <- rnorm(n = nrow(dairy_queen), mean = dqmean, sd = dqsd)

最初の引数は、生成したい数値の数を示しており、nrow()関数を使用してdairy_queenデータセットのメニュー項目の数と同じ数を指定します。 最後の2つの引数は、シミュレートされたサンプルが生成される正規分布の平均と標準偏差を決定します。 シミュレートされたデータセットsim_normの形状と、その正規確率プロットを見ることができます。

  1. sim_normの正規確率プロットを作成してください。すべての点が線上にありますか?このプロットは実データの確率プロットとどのように比較されますか?(sim_normはデータフレームではないため、sample引数に直接入れることができ、data引数は省略できます。)

オリジナルのプロットを正規分布から生成された単一のプロットと比較するよりも優れているのは、以下の関数を使用してより多くのプロットと比較することです。 これは、元のデータに対応するQ-Qプロットを左上隅に、8つの異なるシミュレートされた正規データのQ-Qプロットを表示します。 プロットウィンドウのズームボタンをクリックすると役立つ場合があります。

qqnormsim(sample = cal_fat, data = dairy_queen)
## Warning: `aes_string()` was deprecated in ggplot2 3.0.0.
## ℹ Please use tidy evaluation idioms with `aes()`.
## ℹ See also `vignette("ggplot2-in-packages")` for more information.
## ℹ The deprecated feature was likely used in the openintro package.
##   Please report the issue at
##   <https://github.com/OpenIntroStat/openintro/issues>.
## This warning is displayed once every 8 hours.
## Call `lifecycle::last_lifecycle_warnings()` to see where this warning was
## generated.
  1. 脂肪からのカロリーの正規確率プロットは、シミュレートされたデータに対して作成されたプロットと類似していますか? つまり、プロットは脂肪からのカロリーがほぼ正規分布しているという証拠を提供しますか?

  2. 同じ手法を使用して、マクドナルドのメニューからのカロリーが正規分布から来ているかどうかを判断してください。

正規確率

さて、これで変数が正規分布しているかどうかを判断するための多くのツールが手に入りました。 なぜ気にする必要があるのでしょうか?

統計学者は正規分布について多くのことを知っていることがわかります。 確率変数がほぼ正規分布していると判断すると、その変数に関するあらゆる種類の確率関連の質問に答えることができます。 たとえば、「ランダムに選択されたデイリークイーン製品が脂肪から600カロリー以上を占める確率はどのくらいですか?」という質問を考えてみましょう。

デイリークイーンのメニューからの脂肪からのカロリーが正規分布していると仮定すると(非常に良い近似でも構いません)、Zスコアを計算し、Zテーブル(正規確率テーブルとも呼ばれる)を参照することでこの確率を見つけることができます。 Rでは、これはpnorm()関数を使用して1つのステップで行われます。

1 - pnorm(q = 600, mean = dqmean, sd = dqsd)

pnorm()関数は、与えられた平均と標準偏差を持つ、与えられた値q以下の正規曲線の下の面積を与えることに注意してください。 デイリークイーンのアイテムが脂肪から600カロリー以上を占める確率に関心があるので、その確率から1を引く必要があります。

正規分布を仮定することで、理論的な確率を計算することができました。 確率を経験的に計算したい場合は、600を超える観測値がいくつあるかを特定し、この数を全サンプルサイズで割るだけです。

dairy_queen |>
  filter(cal_fat > 600) |>
  summarise(percent = n() / nrow(dairy_queen))

確率は完全に同じではありませんが、合理的に近いです。 分布が正規分布に近づくほど、理論的な確率はより正確になります。

  1. このデータセットのいずれかのレストランについて、あなたが答えたい2つの確率に関する質問を作成してください。理論的な正規分布と経験的な分布の両方を使用してこれらの確率を計算してください(合計で4つの確率)。どちらの方法の方が一致度が高かったですか?

さらなる演習

  1. 次に、データセットの他の変数について考えてみましょう。 すべてのレストランの中で、ナトリウムの分布が正規分布に最も近いのはどれですか?

  2. ナトリウム分布の正規確率プロットの一部には、段階的なパターンがあることに注意してください。 なぜそうなると思いますか?

  3. ご覧のとおり、正規確率プロットは正規性を評価し、歪度を視覚化するために使用できます。 選択したレストランの総炭水化物量の正規確率プロットを作成してください。 この正規確率プロットに基づいて、この変数は左に歪んでいるか、対称か、右に歪んでいるか? ヒストグラムを使用して、あなたの発見を確認してください。


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