統計的推測

確率密度関数 \( f(x) \) がわからないと、確率を算出することはできません。血糖値の真の分布（関数）は、血糖値の測定を繰り返し行うことによって特定されます。ところが、測定を無限に繰り返すことはできません。そこで、モデルを考えます。過去の経験や現象に基づいて、データのバラツキを再現し、なおかつ、数学的に扱いやすい関数を用います。今、このモデルを表す関数を \( \phi (x) \) とします。

血糖値の具体的な測定値を97.6 mg/dLのように小数点第1位まで表示するとき、小数点第2位以下は四捨五入しています。仮に、無限の精度を持つ測定器を用いれば、小数点以下を限りなく細かく測定できます。このような場合、確率変数は連続的な値をとり、連続型の確率変数といいます。

確率密度関数 \( f(x) \) としてよく使われるものに正規分布があります。これは連続型の確率変数に関する確率分布で、平均値の付近にデータが集まるような分布を表します。 \begin{equation} \phi = \frac{1}{\sqrt{(2\pi \sigma^2)}} \exp \left \{ -\frac{(x-\mu)^2}{2\sigma^2} \right \} \end{equation} \( \phi \) は下図のように、左右対称で釣り鐘型の形状をしています。対称軸は \( \mu \) で正規分布の平均、分布の広がりは \( \sigma \) で正規分布の標準偏差を表します。このような分布は \( N(\mu, \sigma^2) \) と書き表します。平均 \( \mu \) と標準偏差 \( \sigma \) が決まりますと、この正規分布は確定します。

前に出てきた血糖値の表のデータに基づけば、血糖値の平均は96.6 mg/dL, 標準偏差は10.6 mg/dLです。下図の曲線は、\( \mu = 96.6, \sigma = 10.6 \)の正規分布を表しています。ところで、この図を真の分布 \( f(x) \) として、正規分布のモデル \( N(96.6, 10.6^2) \) で近似するとすれば、「高血糖の疑いがある（110 mg/dL 以上）」と診断される確率は約10.2 %です。この場合、10回に1回くらいは「高血糖の疑いがある」と診断されてしまうということです。血糖値の平均は96.6 mg/dL, 標準偏差は10.6 mg/dL というのは本来、標本の標準偏差でしたから、真の分布とは違うという点に注意が必要です。

データにはバラツキがあり、そのバラツキを確率分布としてとらえています。このように、ある事象が起きる確率を、確率分布に基づいて評価することを、統計的推測といいます。

データのバラツキと確率

「データのバラツキは確率的に起きる」と考え、データの背後に隠れている真実を推定します。下図は、先の表で示した、ある人の血糖値を100回測定し、10 mg/dLの階級に分けて描いた頻度グラフです。例えば、血糖値が80～90 mg/dLｎ範囲で測定された回数を、棒の長さでプロットします。このようなグラフはヒストグラムと呼ばれます。全体の測定回数100で割ると、棒グラフの面積は1に正規化されます。

血糖値の測定回数を増やして、階級幅を狭めると、ヒストグラムは図の曲線のような分布に近づきます。この曲線は正の値をとり、横軸との間で囲まれた領域の面積は1となります。つまり、 \begin{equation} f(x) \geq 0 , \int_{-\infty}^{+\infty} f(x) dx = 1 \end{equation} です。これを満たす関数 \( f(x) \) は確率密度関数と呼ばれます。この人の血糖値はバラついていますが、そのバラツキは確率密度関数に従っているということです。この血糖値のデータを変数 \( X \) とすれば、確率密度関数に従って分布する確率変数ということができます。

血糖値が区間 \( (a, b) \) に入る確率は、\( (a, b) \) の範囲の頻度を足せばよいです。確率密度関数で言えば、その範囲の領域の面積を求めればよいのです。この場合、区間 \( (a, b) \) と横軸、\( f(x) \) で囲まれた領域の面積を求めればよいということです。これは、確率変数 \( X \) が区間 \( (a, b) \) に入る確率となります。 \begin{equation} P(a < X < b) = \int_{a}^{b} f(x) dx \end{equation}

これは以前も書いたのですが、確率変数というのはわかりにくい概念です。実際に試行が実施されるまでどのような値が得られるかわからない変数が、確率変数です。試行が実施されて値が得られたら、その値は確率変数の実現値です。確率分布は、確率変数がどのような値を取りやすく、どのような値を取りにくいかを示しています。言い換えれば、確率変数がどのような値をどんな確率で取るかを表すということです。確率変数は、試行を行ってはじめて値（実現値）が決まります。その値に再現性があるわけではなく、同じ試行を繰り返しても同じ値が得られるとは限りません。

統計的な処理

最近、血糖値が上がってきました。もうすぐ健康診断があるんですが、血糖値が基準値を上回るのではないかと心配で、糖質制限と運動を進めています。一般には、高血糖には食事と運動と言われますが、その効果はいかほどなのでしょうか。

例えば、毎夕、食事の後にジョギングを始めたとします。そうしたら、2, 3か月のうちに血糖値が下がってきました（とします）。この血糖値の低下は、運動（ジョギング）の効果で下がったと言えるでしょうか。

「毎日、走ってるんだろう？！血糖値が下がったのなら、「運動の効果」が出たのに決まっている！」

一般には、このように考える人が多いのではないでしょうか。ところが、これは本当に「運動の効果」なのでしょうか。と言うのも、測定にはバラツキがあるはずで、「運動の効果」がなくてもバラツキがあるために、そのときたまたま血糖値が下がったのかも知れません。「運動の効果」によって「血糖値が下がった」というのは早計ではないでしょうか。

このような場合、統計処理を行います。「血糖値が低下したのは運動の効果であり、測定のバラツキによるものではない」ということを知る手段としてp値があります。これは、「もし、運動の効果（処置効果）がないとしたとき、測定のバラツキによって血糖値（評価指標）が、測定された値より小さな値をとる確率」として考えます。

この確率が5%（20回のうち1回）より小さくなるような値が出るならば、バラツキだけの要因で血糖値がこのような小さな値をとる可能性は小さいはずです。したがって、運動の効果があったと推論します。つまり、

• p値 \( < 0.05 \) のとき、統計的に処置効果がある
• p値 \( \geq 0.05 \) のとき、統計的には処置効果があるとは言えない

ここで、「統計的に処置効果がある」というのは統計的な扱いだという点に注意します。

何か測定したとき、データが得られます。繰り返し同じ測定をしたとき、同じ測定値が得られるとは限りません。以下の表は、同じ人の血糖値を測定した結果です。空腹時血糖値を毎日、10日間測定したものです。同一人物、同一条件で測定したにもかかわらず、測定値はバラついています。

表：血糖値の測定（単位：mg/dL）

日にち	0	1	2	3	4	5	6	7	8	9
血糖値	87	96	89	104	88	121	102	88	99	92

血糖値は血液中のブドウ糖の濃度で表され、人によって異なります。しかし、上の表は同一人物でも、血糖値は毎日変動しています。この人の真の血糖値はあるでしょうが、測定値はその真の値に何らかの誤差が加わったものと考えられます。これが、データのバラツキとして出現します。誤差の要因には様々なものが含まれ、一般には特定するのが難しいです。

ところが、健康診断では1回しか測定しないことが多く、その測定で121 mg/dLと出たら「高血糖の疑いあり」と診断されるでしょう（空腹時血糖値110 mg/dL以上で異常）。このように測定データにはバラツキがありますから、精密検査を受けて診断を確定するプロセスがあるのです。したがって、単にデータに基づいて判断するというだけではなく、そのバラツキを認識してデータの背後にある真実を見抜かないといけません。統計学はデータから効率的に情報を抽出し、真実を精度よく推定する方法を提供します。

特異値分解

特異値分解についてお話しします。固有値分解が対象としているのは正方行列でした。特異値分解では、これを一般の行列に拡張します。つまり、行列 A に対して、直交行列 U, V、対角行列 S を用いて \begin{equation} \rm A = U S V^t \end{equation} のように分解します。 このとき、 \begin{equation} \rm S = \begin{pmatrix} s_1 & & \\ & s_2 & \\ & & \ddots \\ \end{pmatrix}, \\ \rm U = \begin{pmatrix} \vec{u_1} & \vec{u_2} & \cdots & \\ \end{pmatrix} \\ \rm V = \begin{pmatrix} \vec{v_1} & \vec{v_2} & \cdots & \\ \end{pmatrix} \end{equation} です。\( \rm S \) は特異値と言います。分解するときに使う直交行列はU, Vの2個になりますが、「正方行列の」という制限が外れます。従って、どんな行列にも使えます。ここで注意すべき点は、対角行列 S の行と列の数が異なるということです。対角からはみ出る要素はゼロになります。また、直交行列 U, V の要素数は異なります。

行列 A は正方行列である必要はありませんが、\( \rm AA^t \) は常に正方行列（対称行列）になります。ということは、この正方行列 \( \rm AA^t \) に対して固有値分解できるということです。つまり、対角行列 D として \begin{equation} \rm AA^t = U D U^t \end{equation} と変形できます。一方で、 \begin{equation} \rm AA^t = U S V^t (U S V^t)^t = U S V^t V S^t U^t = U SS^t U^t \end{equation} ですから、\( \rm D = SS^t \)です。同様に、 \begin{equation} \rm A^t A = V S^t S V^t \end{equation} です。

\begin{equation} \rm S^t S = \begin{pmatrix} s_1^2 & & \\ & s_2^2 & \\ & & \ddots \\ \end{pmatrix} \end{equation} ですから、\( \rm A^t A \) を求めて、その固有値を計算し、その正の平方根をとれば、特異値が求まるということです。

固有値分解

固有値と固有ベクトルの話は少し前にやりました。今回は固有値分解についてです。固有値分解は実質的に対角化そのものです。で、対角化とは何だったかと言いますと、正方行列 A に対して、直交行列 P, 対角行列 D として \begin{equation} \rm {A = P D P^t} \end{equation} のような形にすることです。このとき、 \begin{equation} \rm D = \begin{pmatrix} \lambda_1 & & \\ & \lambda_2 & \\ & & \ddots \\ \end{pmatrix}, \\ \rm P = \begin{pmatrix} \vec{v_1} & \vec{v_2} & \cdots & \\ \end{pmatrix} \end{equation} です。ただし、\( \lambda \) は固有値、\( \vec{v} \) は固有ベクトルです。D は固有値が対角成分に並んだ形になっていますから、対角化です。

で、 \( \rm {A = P D P^t} \) という形に変換したことは、正方行列 A を P の要素（ P に含まれる直交ベクトル）それぞれの和で表したことになります。つまり、 \begin{equation} {\rm A} = \lambda_1 \vec{v_1} \vec{v_1^t} + \lambda_2 \vec{v_2} \vec{v_2^t} + \cdots \end{equation} です。固有値は各項の重みを表していることになりますから、ゼロに近い固有値は無視して計算してしまっても、ある程度の精度が保てます。これは次元削減したことになり、計算コストを下げたりデータ量を減らしたりできます。あとは、\( \rm A \)の累乗の計算が簡単になります。

Google Colaboratoryを使えるようにする

Google Colaboratoryを使おうと思って、net検索すると、Google Driveから「アプリを追加」で対応するようなことが記載されています。で、その通りやって、Google Colaboratoryを登録するのだけど、なぜか関連のフォルダもファイルもできないし、一向に使えるようになりません。

次のようにすると、自動的にGoogle Driveにフォルダやファイルが作成されて、Google Colaboratoryが使える状態になりました。

1. chromeでGoogleにログインする（IDとPWを入力する）
2. chromeから直接Google Colaboratoryのサイトにアクセスする
3. メニューの「ファイル」タブをクリックし、「Python3の新しいノートブック」をクリックする

これで、Google Driveに「Colab Notebooks」というフォルダができ、その中に「Untitled0.ipynb」というファイルができているはずです。

MinGWでOpenCVを使えるようにする

前回、MSYS2をインストールして、gcc関連のパッケージを導入しましたので、C++のプログラムがコンパイルできるはずです。で、以下のようなC++プログラムを「test01.cpp」というファイル名で用意したとします。

#include <iostream>
using namespace std;
int main(){
  cout << "Hello, World." << endl;
  return 0;
}

MinGW 64bitのshellを起動し、test01.cppを作ったディレクトリに移動します。

cd f:\OpecCV

で、以下のコマンドを入力します。

gcc test01.cpp -o test01 -lstdc++

そうしますと、「test01.exe」という実行ファイルができます。これを実行しますと、「Hello, World」が表示されます。

次に、OpenCVを導入します。OpenCVは、コンピュータ・ビジョンに関する機能を持つライブラリです。MSYS2ではOpenCVパッケージが用意されていて、以下のようにするとインストールできます。

pacman -S mingw-w64-x86_64-opencv

で、以下のようなC++プログラムを「test02.cpp」というファイル名で用意したとします。

#include <opencv2/core/core.hpp>
#include <opencv2/highgui.hpp>
#include <iostream>
using namespace cv;
using namespace std;
int main()
{
    Mat img = cv::imread("image.png");
    imshow("image", img);
    waitKey(0);
    return 0;
}

このプログラムをgccでコンパイルするのですが、設定が色々面倒ですので、以下のようなMakefileを作ります。

#Makefile
INC = -I "C:\msys64\mingw64\include\opencv4"
DIR =  -L "C:\msys64\mingw64\bin"
FLAGS = -lstdc++ -lopencv_core410 -lopencv_highgui410 -lopencv_imgcodecs410
OBJS = test02

main:
 gcc $(INC) $(DIR) test02.cpp -o $(OBJS) $(FLAGS)

で、以下のようにmakeすると、実行ファイルができます。

make

この実行ファイルを実行すれば、画像が表示されます。