サーバ構築との出会い

 現在、学部1年生対象の「オペレーティングシステム」という授業を担当しています。この授業では、クラウド環境(AWS)を活用して、サーバ構築を学んでいく授業です。僕自身、サーバ構築を学んだというより、色々なことを頼まれたり、やったりしながら経験を積んでいきました。

(1) 最初の出会い：Linuxと研究室サーバ

僕が学生だった頃(2000年前後)はLinuxのサーバ利用が広がっていた時代でした。大学にUnixが整備されていましたが、ワークステーション室にあるので、大学にろくに行っていなかった僕にとっては遠い存在でした。その中で、Linuxと出会い、無料で使えるサーバに適したOSということで、少し興味を持っていました。その頃の大学のネットワークは、全部IPアドレスはグローバルアドレスを使うなど、今のネットワーク環境では考えられないものでした。僕は、研究室のネットワークを管理しようと、プリンタサーバ、ファイルサーバ(SAMBA)、DHCPサーバなど研究室内で必要なサーバをLinuxで構築しました。これが僕のサーバ構築の最初の経験となりました。

(2) 自宅でのサーバ構築

大学院から社会に出て働き出した頃、自分でHPを作成してみようと思い、せっかくなのでサーバも立ててみようと思って、固定IPやドメインを契約しました。当時のマウスコンピュータで、サーバとするPCも購入し、自宅でサーバを立てて運用していました。今思えば、この１から本番用のサーバを立てて、身近な場所で運用していたことは、良い経験になりました。ただし、自宅でのサーバ運用は、電気代、騒音など苦労も多かったです。

(3) クラウド・コンピューティングとの出会い

2010年前後からクラウド・コンピューティングが急激に普及して、今までのオンプレミスの環境からクラウド環境への移行が進みました。僕自身、静岡大学情報基盤センターで客員教員として研究に関わることができ、当時のセンター長の井上春樹先生が積極的にクラウド活用を進めていたこともあり、色々な研究に関わらせてもらいました。このクラウド環境の出現で、様々な情報基盤の効率化や最適化が図られたと思います。このクラウド環境では、「いつでも」「どこでも」「すぐに」サーバを構築することができ、AWSのEC2は研究だけでなく、授業でも活用しました。学会のHPも今までオンプレミスだったものを、クラウドに移行したのを覚えています。今では当たり前のクラウド環境ですが、この出会いは僕の研究の幅を大きく広げてくれました。

(4) 情報基盤系センター長としての経験

僕自身、8年間情報系基盤センターの副センター長・センター長を担当し、大学のネットワークやサーバ環境を大規模に構築しました。大学の基幹サーバをデータセンタに移行し、ネットワークもデータ活用面まで考えた設計にしました。この取り組みは学内外で評価して頂き、色々なところでお話させてもらいました。ここで内容を述べると大変長くなってしまうのでやめておきますが、いつか書き残しておきたいと思います。この情報基盤構築を経験したことで、現在やっているデータサイエンスの環境を構築することができています。強固な情報基盤を構築し、日々の運用の中で発生するデータを築することができ、さらにそのデータを活用し、情報基盤にフィードバックすることができる環境はデータサイエンスの基盤となります。

(5) 今後の取り組み

僕自身の経験を考えてみると、自宅で1からサーバを立てて運用した経験がベースになり、クラウド環境の活用を行い、大学の情報基盤構築といった大規模な環境まで繋げることができました。クラウド環境が当たり前の現在では、中々サーバというものが目にすることもなく、身近なものではありませんが、物理的なモノを見ることでイメージも湧きやすいと思います。機会があれば、raspberry piなどのマイコンで十分ですので、手元でのサーバ構築を行い、目にみえる経験を得てください。その上でクラウド環境を使うと、サーバ構築の見え方が違うかもしれません。

サーバ環境も日々技術進歩があります。今の知識・技術に固執することなく、生成AIも取り込みながら、次のサーバ環境を考えていきたいと思います。

ベイジアンネットワークの基本1 : 簡単な例

ベイジアンネットワークの簡単な例を通して、どのような仕組みかを確認していきます。

1. 例題1
藤田一弥. (2015). 見えないものをさぐる―それがベイズ ツールによる実践ベイズ統計. 株式会社 オーム社. P99 4-2 例題2を利用して、例題を見ていきます。以下のネットワークを考えます。W(Weather), R(Rain), S(Sprinkler), G(Ground)となっていて、天気の良し悪しで、地面が濡れているかどうかを、雨またはスプリンクラーの影響を入れてネットワークができています。

この例題を通して実施したいことは、以下の3つです。

(1) ベイジアンネットワークのモデル(構造)を作る(ネットワークの形、尤度)

(2) 前方への確率(条件付き確率)を求める (P(G|W=w))など、グラフの左から右への確率

(3) 後方への確率(事後確率)※これがメインのもの、グラフの右から左への確率

今回は、pgmpyを使ってベイジアンネットワークでの分析をするので、インストールをしてください。

!pip install pgmpy

(1) ベイジアンネットワークのモデル(構造)を作る(ネットワークの形、尤度)

まず、例題で与えられたネットワークを作成していきます。また、推移確率も追加していきます。この推移確率を作成するときが間違えやすく、行がvariable(確率変数)、列がevidence(条件)の組み合わせになっているので、推移確率が今回の例題のように与えられた場合は、注意をしてください。今回の推移確率とWの事前分布はソースコードを確認してください。例えば、[[P(S=0|W=0), P(S=0|W=1)], [P(S=1|W=0), P(S=1|W=1)]]のように与えます。

(2) 前方への確率(条件付き確率)を求める (P(G|W=w))など、グラフの左から右への確率

今回は、「P(G = 1| W = 1)を計算：天気が晴れという条件のもとで、芝が濡れている確率」を求めていきます。つまり天気(W=1)が与えられたときのGの分布を求めていきます。

(3) 後方への確率(事後確率)※これがメインのもの、グラフの右から左への確率

今回は、「確率的推論(事後確率) P(W = 1 | G = 1) 芝が濡れているとき、晴れである確率」を求めていきます。「P(W = 1 | G = 1): 0.42538461538461536」が得られます。

2. 例題2 : 涌井良幸, & 涌井貞美. (2012). 史上最強図解 これならわかる! ベイズ統計学. ナツメ社. P170 例題2を参照

もう一つ例題を見ていきます。

この図のようなネットワークを作成します。

このとき、B(Burglar)、E(EarthQuake)、A(Alert)、P(Police)、S(Security)といったノードになっています。泥棒(B)が入ったとき、アラート(A)がなり、セキュリティ会社(S)に連絡がいくという流れです。

(1) ベイジアンネットワークのモデル(構造)を作る(ネットワークの形、尤度)

これは例題で与えられているので、尤度(推移確率)の構造に注意して作成していきます。

(2)前方への確率(条件付き確率)を求める (P(G|W=w))など、グラフの左から右への確率

これは演習としますので、やってみてください。

(3) 後方への確率(事後確率)※これがメインのもの、グラフの右から左への確率

「P(B | S) セキュリティ会社に連絡がいったとき、泥棒が入った確率」を求めてみます。

P(B = 1 | S = 1): 0.042763793338810904が得られます。

ここまでのソースコード

Bakery Sales Dataset : (2) 日毎の販売回数をポアソン分布に従うとして、パラメタの事後分布を推定

 今回は、日毎の販売回数がポアソン分布に従うとして事後分布を推定していきます。

(1) 共役自然分布を使ったポアソン分布の事後分布推定の概要

ポアソン分布の尤度関数は次のように与えられる：

$$ p(X|\lambda) = \prod_{i=1}^{n} \frac{e^{-\lambda} \lambda^{x_i}}{x_i!} $$

ここで、$ X=(x_1,\dots,x_n) $ は観測データである。

ガンマ事前分布は次のように与えられる：

$$ p(\lambda|a,b) = \frac{b^a}{\Gamma(a)} \lambda^{a-1} e^{-b\lambda} $$

ベイズの定理により、

$$ p(\lambda|X) = \frac{p(\lambda,X)}{p(X)} = \frac{p(X|\lambda)p(\lambda|a,b)}{p(X)} \propto p(X|\lambda)p(\lambda|a,b) $$

$$ = \prod_{i=1}^{n} \frac{e^{-\lambda} \lambda^{x_i}}{x_i!} \cdot \frac{b^a}{\Gamma(a)} \lambda^{a-1} e^{-b\lambda} $$

$$ \propto \lambda^{\left(\sum_{i=1}^{n} x_i + a - 1\right)} e^{-(n+b)\lambda} $$

以上より、事後分布はガンマ分布となり、そのパラメータは次のようになる：

$$ a_{\text{post}} = \sum_{i=1}^{n} x_i + a $$

$$ b_{\text{post}} = n + b $$

したがって、$ \lambda $ の事後分布は次のようになる：

$$ p(\lambda|X) = \frac{b_{\text{post}}^{a_{\text{post}}}}{\Gamma(a_{\text{post}})} \lambda^{a_{\text{post}}-1} e^{-b_{\text{post}}\lambda} $$

つまり、ポアソン分布のパラメタ$\lambda$の事前分布をガンマ分布とおくと、$\lambda$の事後分布は同じくガンマ分布になる。

(2) 販売回数がポアソン分布に従うときのパラメタ推定

ここでは、1日ごとの販売回数を集計し、そのデータを用いてポアソン分布のパラメータのベイズ推定を行います。その後、事後分布を求めてみましょう。まずは、日付ごとの販売回数を集計します。

ここでは、事前分布としてガンマ分布を使用し、事後分布を計算してみましょう。ガンマ分布は、形状パラメータ（$\alpha$）と尺度パラメータ（$\beta$）を持ちます。これらのパラメータを適切に選ぶことが重要です。仮に、$\alpha=\beta=1$（一様事前分布に近い）として計算してみます。

ポアソン分布のパラメータ（$\lambda$）の事後分布をプロットしました。このグラフは、観測データとガンマ分布を事前分布とした場合の$\lambda$の事後確率密度を示しています。

この事後分布から、平均販売回数の最も確からしい値や、その不確実性を評価することができます。例えば、事後分布のピークは$\lambda$の最も可能性の高い値を示し、分布の幅は不確実性の程度を示します。

(3) 得られた$\lambda$の事後分布からのMAP（Maximum A Posteriori）推定
AP推定は事後分布を最大化する$\lambda$の値を見つけることです。ガンマ分布の場合、$a > 1$ のとき、$\lambda$のMAP推定値は次のようになります： 

$$ \hat{\lambda}_{\text{MAP}} = \arg\max_\lambda (p(X|\lambda) \cdot p(\lambda|a,b)) $$ 

$$ \hat{\lambda}_{\text{MAP}} = \frac{a_{\text{post}} - 1}{b_{\text{post}}} $$ 

ガンマ分布の最大値はそのモード（最頻値）となる性質があります。

実際求めてみます。

ポアソン分布のパラメータ（λ）の最大事後確率（MAP）推定値は約59.16です。これは、観測データをもとに、平均販売回数が一日あたり約59回であると推定されることを意味します。この値は、データと事前分布から導き出された最も確からしい販売回数の推定値です。

今度は、$\alpha=2, \beta=1$でやって、先ほどの結果と比較してみます。新しい事前分布に基づいた$\lambda$の最大事後確率（MAP）推定値は約59.16で同じでした。これは、事前分布を変更しても、観測データの影響が大きいため、MAP推定値に大きな変化が見られないことを示しています。この結果は、データが豊富であればあるほど、事前分布の影響が小さくなるというベイズ推定の特性を反映しています。

(4)事後予測分布

事後予測分布（Posterior Predictive Distribution）は、ベイズ統計学の枠組みで使用される重要な概念の1つで、ベイズ推論を行った結果から未知のデータポイントの分布を予測するための確率分布です。事後予測分布は以下のように表されます: 

$$ P(X_{\text{new}}|X) = \int P(X_{\text{new}}|\lambda)P(\lambda|X)d\lambda $$ 

ここで、ポアソン分布のパラメータ $\lambda$ がガンマ分布に従う場合、事後予測分布は負の二項分布（またはポリヤ分布）になります。 

負の二項分布の確率質量関数（PMF）は以下のように表されます: 

$$ P(x_{\text{new}}|a_{\text{post}}, b_{\text{post}}, X) = \frac{\Gamma(x_{\text{new}} + a_{\text{post}})}{x_{\text{new}}!\Gamma(a_{\text{post}})} \left(\frac{b_{\text{post}}}{b_{\text{post}} + 1}\right)^{a_{\text{post}}} \left(\frac{1}{b_{\text{post}} + 1}\right)^{x_{\text{new}}} $$

これをやろうとすると、先ほど求めたパラメタ$\alpha_post=9466, \beta=160 (\alpha=2, \beta=1)$となってしまい、累乗の計算が困難です。以下のように実施します。

最後にMCMCで事後予測分布を可視化します。

ここまでのソースコード

糖尿病：Diabetes Health Indicators Dataset : 02 主成分分析からクラスタリングまで

 糖尿病：Diabetes Health Indicators Dataset : 基本分析01の続きを行なっていきます。今回は主成分分析と主成分分析の結果を用いたクラスタリングです。

1. 主成分分析の実施

まず主成分分析を行なっていきますが、今回のように要素が多い場合、主成分分析で次元削減を行い、新たな軸で分析することは有効です。今回は累積寄与率80%で分析を行います。主成分分析を行い、各主成分の項目ごとのウエイトをデータフレームにして算出します。今回は目的変数に相当する「Diabetes_binary」は除いた21項目で実施しました。80%の累積寄与率でやると14主成分となっています。

これだとはっきりしないので、ヒートマップで色分けしてみます。

演習：このヒートマップをみて、各主成分(PC0 ~ PC4)に名前をつけましょう。

• PC0: この主成分は「GenHlth」（一般的な健康状態）や「PhysHlth」（身体的健康）などの特徴量に高い負荷量を持っています。したがって、この主成分は「全体的健康状態」を反映していると考えることができます。
• PC1: 「Age」（年齢）に高い負荷量を持っており、「Education」（教育）や「Income」（収入）にも影響を与えています。この主成分は「年齢と社会経済的地位」を反映していると言えるでしょう。
• PC2: 「Fruits」（果物の摂取）や「Veggies」（野菜の摂取）に高い負荷量を持っています。この主成分は「食生活の健康性」を表している可能性があります。
• PC3: この主成分は「BMI」（体格指数）に高い負荷量を持っており、「CholCheck」（コレステロール検査）や「GenHlth」（一般的な健康状態）にも影響を与えています。この主成分は「身体的健康とライフスタイル」を反映していると考えられます。
• PC4: 「Smoker」（喫煙）に高い負荷量を持っていることが見られます。また、他の健康関連の特徴量にも影響を与えています。この主成分は「喫煙と関連する健康リスク」を表している可能性があります。

累積寄与率のグラフを書いてみます。

次に主成分0と1を使って、biplotで表示します。まずは「Diabetes_binary」が1のときです。

次に「Diabetes_binary」が0のときです。

PC2も加えて3次元で表示してみます。

演習：biplotの結果からどのようなことが言えるでしょうか？

2. 主成分分析の結果からクラスタリング実施

次にクラスタリングを実施してみます。元データからではなく、主成分分析で求まる射影行列を使ってのクラスタリングを行います。今回はK-means法で5個のクラスタに分類してみます。

この結果をヒートマップで表示してみやすくしてみます。値を整えてから出力します。

演習：それぞれのクラスタの特徴を書き出してみましょう。

• クラスタ0 (「低リスク・健康志向」): このクラスタは全体的に低い糖尿病率と低い健康リスク指標（BMI、高血圧など）を示しており、身体活動が高いことが特徴です。

• クラスタ1 (「中高年・リスクあり」): 高い年齢層に属し、高血圧や高コレステロールの割合が高い。中程度の糖尿病率を持つ。

• クラスタ2 (「高リスク・健康課題」): 最も高い糖尿病率を持ち、BMIや心臓病のリスクも高い。身体活動が低く、健康上の課題が多い。

• クラスタ3 (「アクセス限定・中リスク」): 一般的な健康状態は中程度で、医療へのアクセスが限られている（NoDocbcCostが高い）ことが特徴。中程度の糖尿病率。

• クラスタ4 (「若年・喫煙者」): 比較的若い年齢層で喫煙率が高い。糖尿病の割合は比較的低いが、長期的な健康リスクを持つ可能性がある。

このような方法で、基本的な分析からクラスタリングをすることで、データの内容がわかりやすくなったと思います。

ここまでのコード：https://github.com/smzn/diabetes_binary_health_indicators_BRFSS2015/blob/main/diabetes_binary_health_indicators_BRFSS2015_02.ipynb