BioinfoにおけるPandas,Matplotlibの基礎

この記事は、Pythonで実践　生命科学データの機械学習 (https://www.yodosha.co.jp/yodobook/book/9784758122634/) の内容を含んでいます。

私は現在バイオインフォマティクス研究室に所属する学生です。
勉強した事をアウトプットする場として用いていますため、何卒ご理解のほどよろしくお願いします
(＞人＜;)

マイクロアレイデータをPyTorchを用いで分類

マクロアレイは、DNAやRNA、タンパク質などの分子を網羅的に解析するための高密度な配列技術である。特定の遺伝子の発現量や変異、コピー数変化を測定する際に利用される。

本日は、深層学習ライブラリPyTorchを用いて機械学習分類器の構築を学んだ。その際に、PyTorchの基本的な文法や知識を学んだので記す。

PyTorchチュートリアル

PyTorchでは、計算は「テンソル」というデータ型を用いて行われる。テンソルはスカラー（階数0）、ベクトル（階数1）、行列（階数2）の一般化された概念であり、さらに高次元のデータも表現できる。例えば、フルHD画像（1920×1080ピクセル）のRGB情報を持つ場合、これは階数3のテンソル（1920×1080×3）で表現される。

また、PyTorchはNumpyに似たテンソル演算を提供し、大きな特徴としてGPUを用いた高速演算が可能である。Google Colaboratory上でGPUを利用するには、ランタイム設定を「ハードウェアアクセラレータ」に変更する必要がある。設定後、GPUが正しく利用できるかを確認するコードも紹介されている。

基本的なテンソル計算

配列の書き方

# Numpyのアレイ
import numpy as np

data = [[1, 2], [3, 4]]
np_data = np.array(data)
print(np_data)

#出力結果
[[1 2]
 [3 4]]


# PyTorchのテンソル
import torch

tensor_data = torch.tensor(data)
print(tensor_data)

#出力結果
tensor([[1, 2],
        [3, 4]])


#Numpyのアレイエをテンソルに変換
tensor_from_np = torch.from_numpy(np_data)
tensor_from_np

#出力結果
tensor([[1, 2],
        [3, 4]])

Numpy同様に値の代入ができる。

# 4*4の行列の成分がすべて1であるテンソル
tensor = torch.ones(4, 4)

# GPU上へtensorを動かす
tensor.to(device)
print(tensor)

#出力結果
tensor([[1., 1., 1., 1.],
        [1., 1., 1., 1.],
        [1., 1., 1., 1.],
        [1., 1., 1., 1.]])


# Numpyと似たような指定ができる
# 値の代入
tensor[:, 1] = 0
tensor[1, :] = 3
print(tensor)

#出力結果
tensor([[1., 0., 1., 1.],
        [3., 3., 3., 3.],
        [1., 0., 1., 1.],
        [1., 0., 1., 1.]])

Numpy同様に値の抽出もできる。

# 行や列の抽出
print('0行目: ', tensor[0])
print('1列目: ', tensor[:, 1])
print('最終列: ', tensor[:, -1])

#出力結果
0行目:  tensor([1., 0., 1., 1.])
1列目:  tensor([0., 3., 0., 0.])
最終列:  tensor([1., 3., 1., 1.])

コピーして横に貼り付けることもできる

# dim=1で横に並べる
t1 = torch.cat([tensor, tensor, tensor], dim=1)
print(t1)

#出力結果
tensor([[1., 0., 1., 1., 1., 0., 1., 1., 1., 0., 1., 1.],
        [3., 3., 3., 3., 3., 3., 3., 3., 3., 3., 3., 3.],
        [1., 0., 1., 1., 1., 0., 1., 1., 1., 0., 1., 1.],
        [1., 0., 1., 1., 1., 0., 1., 1., 1., 0., 1., 1.]])

演算式の一例

# テンソルの計算: ここでは行列式
tensor @ tensor

#出力結果
tensor([[ 3.,  0.,  3.,  3.],
        [18.,  9., 18., 18.],
        [ 3.,  0.,  3.,  3.],
        [ 3.,  0.,  3.,  3.]])


# アダマール積
tensor * tensor

#出力結果
tensor([[1., 0., 1., 1.],
        [9., 9., 9., 9.],
        [1., 0., 1., 1.],
        [1., 0., 1., 1.]])


print(tensor, '\n')#改行の意味で\n
# 'add'ではなく'add_'
tensor.add_(5)
print(tensor)

#出力結果
tensor([[1., 0., 1., 1.],
        [3., 3., 3., 3.],
        [1., 0., 1., 1.],
        [1., 0., 1., 1.]]) 

tensor([[6., 5., 6., 6.],
        [8., 8., 8., 8.],
        [6., 5., 6., 6.],
        [6., 5., 6., 6.]])

テンソルは、数値型変数に変換できる

sum_item = tensor.sum().item()
print(sum_item, type(sum_item))

少し理解に時間がかかったので、詳しく解説していく

このコードの目的は、PyTorch のテンソルの全要素の合計を計算し、その値を Python の組み込みデータ型（float または int）として取得することである。

よくある用途としては、深層学習モデルの損失値やスカラー値を取り出してログに記録する場合に使用されます

loss = torch.tensor(1.234)
print(loss.item())  # 1.234

型の確認と変換.item() を使うことでテンソル型から Python の標準型に変換可能。

深層学習の基本知識

深層学習の基礎知識

1. ニューラルネットワークとは

ニューラルネットワークは、脳の神経細胞（ニューロン）をモデル化した計算システムである。入力データを処理して、分類や予測を行うために利用される。

基本構成

入力層: データを受け取る層。
隠れ層: 入力データを処理・変換する層。
出力層: 最終的な結果を出力する層。

計算グラフ

ニューラルネットワークは「計算グラフ」として表されます。ノード（丸）とエッジ（線）で構成される。

ノード: 計算を行う単位。
エッジ: ノード同士を結び、情報（数値）を流す線。各エッジには「重み」という値が付与され、入力情報の重要度を調整する。

2. 全結合層（Fully Connected Layer）

全結合層は、ある層のすべてのノードが次の層のすべてのノードに接続されている構造を指します。

例えば、入力層が3つのノード ((x_1, x_2, x_3)) を持ち、隠れ層が2つのノード ((y_1, y_2)) を持つ場合、それぞれのノード間で重み ((w_{11}, w_{12}, w_{21}, \dots)) を掛け算して次の層に情報を送ります。

計算式

$$
y_1 = w_{11}x_1 + w_{12}x_2 + w_{13}x_3 + b_1
$$
$$
y_2 = w_{21}x_1 + w_{22}x_2 + w_{23}x_3 + b_2
$$

ここで:

(w_{ij}): 重み
(b_j): バイアス（調整用の定数）

3. 活性化関数

活性化関数は、各層の出力に「非線形性」を加えるために使用します。

主な活性化関数

Sigmoid関数:
$$
f(x) = \frac{1}{1 + e^{-x}}
$$

出力を0～1の範囲に変換。
確率を扱う場合に使用。

ReLU関数（Rectified Linear Unit）:
$$
f(x) = \max(0, x)
$$

負の値を0にし、正の値をそのまま通す。
現在の深層学習で広く使用されている。

Tanh関数:
$$
f(x) = \frac{e^x - e^{-x}}{e^x + e^{-x}}
$$

出力を-1から1の範囲に変換。
Sigmoidよりも学習の収束が早いことがある。

4. 損失関数（Loss Function）

損失関数は、ニューラルネットワークの予測結果と実際の値との「誤差」を数値化します。

主な損失関数

バイナリ交差エントロピー（Binary Cross Entropy）
二値分類問題に適した損失関数で、以下の式で表されます。
$$
L = -\frac{1}{n} \sum_{i=1}^n \left[ t_i \log(p_i) + (1 - t_i) \log(1 - p_i) \right]
$$

ここで:

(t_i): 実際のラベル（0または1）
(p_i): モデルの予測値（確率）

平均二乗誤差（Mean Squared Error; MSE）
回帰問題に使用される損失関数です。
$$
L = \frac{1}{n} \sum_{i=1}^n (y_i - \hat{y}_i)^2
$$

ここで:

(y_i): 実際の値
(\hat{y}_i): 予測値

5. 学習プロセス

ニューラルネットワークの学習は以下の手順で行われます。

(1) フォワードプロパゲーション

入力データをモデルに流し、出力結果を得ます。

(2) 損失計算

損失関数を用いて、予測結果と正解との誤差を計算します。

(3) バックプロパゲーション（誤差逆伝播法）

損失を基にして各重みやバイアスを調整します。このとき、勾配降下法というアルゴリズムを使用します。

(4) 勾配降下法

損失関数の値を最小にするようにパラメータを更新します。重みの更新は以下の式で行われます。
$$
w \leftarrow w - \eta \frac{\partial L}{\partial w}
$$

ここで:

(\eta): 学習率（重みをどれだけ変化させるかを決める値）
(\frac{\partial L}{\partial w}): 損失関数の重みに対する勾配

6. 学習率と収束

学習率: モデルが一度のステップでどれだけ進むかを決めるパラメータ。
- 学習率が小さい場合: 収束が遅い。
- 学習率が大きい場合: 最適な値を飛び越えてしまう。

バイナリー交差エントロピー

バイナリー交差エントロピー（Binary Cross-Entropy, BCE）は、2クラス分類問題でよく使われる損失関数です。
例えば、猫か犬かを判定する画像分類や、正常か異常かを判断する異常検知タスクで使用されます。

この記事では、NumPyを使った自作の計算関数と、PyTorchの組み込み関数を使った計算方法を紹介します。

バイナリー交差エントロピーの数式

バイナリー交差エントロピーの損失関数は以下の数式で表されます：

$$
\text{BCE} = -\frac{1}{N} \sum_{i=1}^N \left[ y_i \log(\hat{y}_i) + (1 - y_i) \log(1 - \hat{y}_i) \right]
$$

\( N \): サンプル数（バッチ内のデータ数）
\( y_i \): 実際のラベル（0 または 1 の値）
\( \hat{y}_i \): モデルの予測確率（\( 0 \leq \hat{y}_i \leq 1 \)）
\( \log \): 自然対数

直感的な理解

予測が正確な場合（損失が小さい）:
- $y_i$ = 1 のとき、$ \hat{y} $ が 1 に近いほど良い
- $y_i$ = 0 のとき、$ \hat{y} $ が 0 に近いほど良い
誤った予測（損失が大きい）:
- 予測確率が実際のラベルと大きく異なるほど損失が増加

PyTorchを使ったバイナリー交差エントロピーの計算

PyTorchには、バイナリー交差エントロピーを簡単に計算できる組み込み関数 torch.nn.BCELoss が用意されている。

import torch
import torch.nn as nn

# サンプルデータ
y_true = torch.tensor([1, 0, 1, 0], dtype=torch.float32)  # 実際のラベル
y_pred = torch.tensor([0.9, 0.1, 0.8, 0.4], dtype=torch.float32)  # モデルの予測確率

# PyTorchの損失関数（Binary Cross-Entropy）
bce_loss = nn.BCELoss()

# バイナリー交差エントロピーを計算
loss = bce_loss(y_pred, y_true)
print("Binary Cross-Entropy Loss (PyTorch):", loss.item())

##実行結果
Binary Cross-Entropy Loss (PyTorch): 0.21616187691688538

y_true: 実際のクラスラベル（正解データ）
[1, 0, 1, 0] → 1番目と3番目のサンプルはクラス1（ポジティブ）、2番目と4番目のサンプルはクラス0（ネガティブ）。
y_pred: モデルが出力した確率値
[0.9, 0.1, 0.8, 0.4]

例えば、最初の値 0.9 は「1番目のデータがクラス1である確率は 90%」という予測を意味する。

bce_loss = nn.BCELoss()

nn.BCELoss() は PyTorch のバイナリー交差エントロピー損失関数である
前提: y_pred はシグモイド関数を適用済みの確率値である必要があります。（0 ≤ y_pred ≤ 1）

loss = bce_loss(y_pred, y_true)
print("Binary Cross-Entropy Loss (PyTorch):", loss.item())

各データについて、次の式を計算し、全体の平均を求める。

$$
-\left( y \log(\hat{y}) + (1 - y) \log(1 - \hat{y}) \right)
$$

ここで、

y は 実際のクラスラベル（0 または 1）
$ \hat{y} $は モデルの予測確率（0 ≤ $ \hat{y} $ ≤ 1）

この値を全データに対して求め、平均することで バイナリー交差エントロピー損失 を計算する。

この損失値の解釈

0.216 という値は比較的小さいため、予測は比較的正確であるといえる。
もし y_pred = [0.5, 0.5, 0.5, 0.5] などのいい加減な確率であれば、損失は 0.693（= log(2)）付近になる。
より低い値（例えば 0.01）になれば、モデルの精度が高いことを意味する。

誤差逆伝播法（Backpropagation）とは？

1. 誤差逆伝播法の概要

誤差逆伝播法（Backpropagation, バックプロパゲーション）は、ニューラルネットワークの学習において、損失（誤差）をネットワークの各層に遡って伝播し、重みを最適化するためのアルゴリズムである。

誤差逆伝播法を使うことで、ニューラルネットワークは「どの重みをどの程度調整すれば損失が減るか」を学習することができる。

2. 誤差逆伝播法の基本的な流れ

誤差逆伝播法の流れは、以下の 順伝播（Forward Propagation）→ 損失計算 → 逆伝播（Backpropagation）→ 重み更新 の4ステップで進む。

① 順伝播（Forward Propagation）

入力データをニューラルネットワークに通し、出力を計算する。

$$
\hat{y} = f(WX + b)
$$

$X$ ：入力データ
$W$ ：重み（学習対象）
$b$：バイアス
$f$ ：活性化関数（例：シグモイド, ReLU）
$\hat{y}$ ：予測値

② 損失関数を計算

出力 $\hat{y}$と正解ラベル $y$を比較し、損失（エラー）を求める。

$$
L = \text{Loss}(\hat{y}, y)
$$

例：二乗誤差（MSE）
$$
L = \frac{1}{2} (\hat{y} - y)^2
$$

例：バイナリー交差エントロピー（BCE）
$$
L = -y \log(\hat{y}) - (1 - y) \log(1 - \hat{y})
$$

③ 誤差逆伝播（Backpropagation）

損失 ( L ) を各層に遡り、各パラメータ（重みやバイアス）の勾配（変化量）を計算する。

偏微分を利用して、重み ( W ) の勾配を求める：

$$
\frac{\partial L}{\partial W}
$$

この勾配を求めるために、連鎖律（Chain Rule） を用いる：

$$
\frac{\partial L}{\partial W} = \frac{\partial L}{\partial \hat{y}} \cdot \frac{\partial \hat{y}}{\partial W}
$$

④ 勾配降下法（Gradient Descent）でパラメータ更新

求めた勾配を使って、勾配降下法（Gradient Descent） により、重み $ W $ を更新する。

$$
W := W - \alpha \frac{\partial L}{\partial W}
$$

$\alpha$学習率, Learning Rate）：どれくらいの速さで更新するかを決めるハイパーパラメータ
勾配が大きい場合、重みを大きく変更し、勾配が小さい場合、微調整する

3. Python（NumPy）で誤差逆伝播法を実装

単純な1層ニューラルネットワークの逆伝播

import numpy as np

# シンプルな単層ニューラルネットワーク
def sigmoid(x):
    return 1 / (1 + np.exp(-x))

def sigmoid_derivative(x):
    return x * (1 - x)  # シグモイド関数の微分

# 訓練データ
X = np.array([[0, 0], [0, 1], [1, 0], [1, 1]])  # 入力
y = np.array([[0], [1], [1], [0]])  # XORの正解ラベル

# 重みとバイアスの初期化
np.random.seed(42)
W = np.random.randn(2, 1)  # 2入力 → 1出力の重み
b = np.random.randn(1)

# 学習率
alpha = 0.1

# 1000回の学習
for epoch in range(1000):
    # 順伝播
    z = np.dot(X, W) + b
    y_pred = sigmoid(z)

    # 損失計算（二乗誤差）
    loss = np.mean((y_pred - y) ** 2)

    # 逆伝播
    dL_dy = 2 * (y_pred - y)  # 損失関数の微分
    dy_dz = sigmoid_derivative(y_pred)  # シグモイドの微分
    dz_dW = X.T  # 入力Xに関する勾配

    # 重みの勾配計算
    dL_dW = np.dot(dz_dW, dL_dy * dy_dz)
    dL_db = np.sum(dL_dy * dy_dz)

    # 勾配降下法による更新
    W -= alpha * dL_dW
    b -= alpha * dL_db

    # 100回ごとに損失を表示
    if epoch % 100 == 0:
        print(f"Epoch {epoch}: Loss = {loss:.4f}")

# 学習後の出力を確認
print("\n学習後の予測:")
print(sigmoid(np.dot(X, W) + b))

PyTorchを使えば、逆伝播（勾配計算）は backward() を呼び出すだけで自動計算される。

import torch
import torch.nn as nn
import torch.optim as optim

# データセット
X = torch.tensor([[0, 0], [0, 1], [1, 0], [1, 1]], dtype=torch.float32)
y = torch.tensor([[0], [1], [1], [0]], dtype=torch.float32)

# モデル定義（単層ニューラルネットワーク）
model = nn.Sequential(
    nn.Linear(2, 1),  # 2入力 → 1出力
    nn.Sigmoid()
)

# 損失関数と最適化アルゴリズム
criterion = nn.MSELoss()  # 二乗誤差
optimizer = optim.SGD(model.parameters(), lr=0.1)  # 勾配降下法（SGD）

# 学習ループ
for epoch in range(1000):
    optimizer.zero_grad()  # 勾配を初期化

    # 順伝播
    y_pred = model(X)

    # 損失計算
    loss = criterion(y_pred, y)

    # 逆伝播（勾配計算）
    loss.backward()

    # 勾配降下法で重み更新
    optimizer.step()

    if epoch % 100 == 0:
        print(f"Epoch {epoch}: Loss = {loss.item():.4f}")

# 学習後の予測
print("\n学習後の予測:")
print(model(X).detach().numpy())

実行結果は以下のようになった。

Epoch 0: Loss = 0.2727
Epoch 100: Loss = 0.2716
Epoch 200: Loss = 0.2706
Epoch 300: Loss = 0.2696
Epoch 400: Loss = 0.2686
Epoch 500: Loss = 0.2676
Epoch 600: Loss = 0.2667
Epoch 700: Loss = 0.2657
Epoch 800: Loss = 0.2648
Epoch 900: Loss = 0.2639

学習後の予測:
[[0.4991]
 [0.5060]
 [0.5072]
 [0.5138]]

結果の解釈

損失が徐々に減少している
Epoch 0: Loss = 0.2727 から Epoch 900: Loss = 0.2639 まで、少しずつ減っている
これはモデルが学習し、誤差（Loss）が改善されていることを示している
学習後の出力がXORの期待値と異なる。

本来 XOR の期待される出力は [0, 1, 1, 0] になるべき
しかし、学習後の y_pred は [[0.4991], [0.5060], [0.5072], [0.5138]] とほぼ 0.5 に近い値になっている
これは単層ニューラルネットワークではXORを学習できないため、非線形性が足りないことが原因である。

改善するには？

隠れ層を追加する
単層のモデルではXORのような非線形問題を解くのは難しい
例えば nn.Linear(2, 3) → nn.ReLU() → nn.Linear(3, 1) のように隠れ層を追加すると学習が進む
活性化関数を変更する
Sigmoid では勾配消失が起きやすいので、ReLU を試すのも良い

ここで隠れ層、活性化関数ReLUを作ってリトライしてみる。

import torch
import torch.nn as nn
import torch.optim as optim

# XORデータ
X = torch.tensor([[0, 0], [0, 1], [1, 0], [1, 1]], dtype=torch.float32)
y = torch.tensor([[0], [1], [1], [0]], dtype=torch.float32)

# 隠れ層を含むニューラルネット
model = nn.Sequential(
    nn.Linear(2, 3),  # 2入力 → 3ノード
    nn.ReLU(),
    nn.Linear(3, 1),
    nn.Sigmoid()
)

# 損失関数と最適化
criterion = nn.MSELoss()
optimizer = optim.SGD(model.parameters(), lr=0.1)

# 学習ループ
for epoch in range(10000):  # エポック数を増やす
    optimizer.zero_grad()
    y_pred = model(X)
    loss = criterion(y_pred, y)
    loss.backward()
    optimizer.step()

    if epoch % 1000 == 0:
        print(f"Epoch {epoch}: Loss = {loss.item():.4f}")

# 学習後の予測
print("\n学習後の予測:")
print(model(X).detach().numpy())

修正後の結果は、このようになった。

Epoch 0: Loss = 0.2674
Epoch 1000: Loss = 0.2502
Epoch 2000: Loss = 0.2067
Epoch 3000: Loss = 0.1392
Epoch 4000: Loss = 0.0803
Epoch 5000: Loss = 0.0431
Epoch 6000: Loss = 0.0236
Epoch 7000: Loss = 0.0135
Epoch 8000: Loss = 0.0081
Epoch 9000: Loss = 0.0049

学習後の予測:
[[0.0134]
 [0.9885]
 [0.9878]
 [0.0122]]

以下のように XORの正解に近い予測になったことがわかる。

活性化関数ReLU

ReLU（Rectified Linear Unit）とは？

ReLU（Rectified Linear Unit）は、ニューラルネットワークでよく使われる 活性化関数 です。
「負の値は0にし、正の値はそのまま通す」 というシンプルな動作をします。

数式

$$
f(x) = \max(0, x)
$$

グラフ

ReLUの関数は以下のように負の部分は0、それ以外は直線になります。

なぜReLUを使うのか？

1. 勾配消失問題を防ぐ

Sigmoid や Tanh は 入力が大きくなると勾配が0に近くなる（学習が進まない）
ReLU は 正の部分では勾配が1のまま なので、学習が止まりにくい！

2. 計算が速い

ReLU(x) = \max(0, x) は 指数関数を使わないので計算コストが低い
Sigmoid や Tanh は指数計算が必要で、計算負荷が高い

ReLUの欠点

1. 「死んだReLU（Dead ReLU）」問題

( x < 0 ) のとき、勾配が0 なので、学習が進まなくなる
これを防ぐために Leaky ReLU（負の値も小さく残すバージョン）がある

2. 出力が大きくなると発散する

Sigmoid は 0~1 の範囲に制限されるが、ReLUは制限なし
重みの初期化（He初期化） や バッチ正規化（Batch Normalization） と組み合わせると安定する