生成AIで画像分類

目から鱗の文献でした。自分の頭の固さを痛感します。
画像キャプショニングは画像そのものよりも多くを語る

画像キャプショニング（image-to-text）の技術を活用して、入力画像に対する説明文を自動生成し、その情報から画像を分類する取り組みでした。これまで、生成AIの発展は自分の仕事には直接関係ないと考えていましたが、その認識は大きな誤りでした。反省です。

Transformerと同様の仕組みが採用されているようで、事前学習モデルも配布されています。ファインチューニングによって精度を出し易くなっているのでしょう。
UAVなどで撮影した動画から構造物の異常を抽出するアプローチにおいては画像解析やセグメンテーションといった従来技術が主流でした。が、生成AIを活用すれば、既に正確な異常検知が可能になっているのではないでしょうか。これまで画像を扱ってきた方々の技術の延長とはやや遠い手法ですので、なかなか思いつかないでしょうね。https://phreeqc.blogspot.com/2025/02/blog-post_26.html
道路分野では、ドライブレコーダーから画像を受けて、リアルタイムで注意を促すなどといった使い方もできそうです。被災後に走るだけで異常個所の状況が写真とともにサーバーに集約される、というようなシステムも組めそうです。自動車メーカーさんなら既に手を付けられているかもしれません。

異常の有無だけなら実装は容易でしょうが、状況をより詳しく説明させるとなると、学習させるテキストデータが重要でしょう。画像の特徴（名称、異常箇所など）を適切に盛り込み、状況を明確に教える必要があります。プロンプトを利用する場合も同様でしょう。さらに、危険度の評価基準（危険度A,B,Cなど）や、対応策のレベル（緊急対策、監視など）を具体的に記述しておくべきでしょうね。そうすることで、生成AIが汎用性を保ちながらも、個々の状況に応じた適切な分析結果を出力できるようになると思います。

適当なデータセットが手元にないので真の実力を把握することはできませんが、単純なマルチクラス分類ならいくつかあります。ひとまずそれらで生成AIの実力を見てみましょうか。

**********************************
20250331 追記
BLIP を利用したのですが、Transformer ベースのアンサンブルモデルの結果を3回のトライで超えました。しかも、学習にそれほど時間がかかりません。生成AI、恐るべし。

勾配消失と回避策

最後に勾配消失。昔はこれがネックだったようです。

1. シグモイド関数の微分の最大値

シグモイド関数の微分は、その定義から以下のように導出される：

σ′(x) = σ(x)(1 − σ(x))

σ(x) の値は常に 0 と 1 の間にある。ここで、関数 g(u) = u(1 − u) （ただし u = σ(x)）が u ∈ (0, 1) で最大になる u の値を求めるため一階微分を計算すると：

g(u) = u(1 − u)   ⇒   g′(u) = 1 − 2u

g′(u) = 0 となる点を求めると：

1 − 2u = 0   ⇒   u = 0.5

この u = 0.5 における g(u) の値を求めると：

g(0.5) = 0.5 × (1 − 0.5) = 0.25

よって、σ'(x) の最大値は 0.25 となる。つまり、シグモイド関数の微分は入力 x に依存して変化するが、その最大値は 0.25 である。なお、σ(x) が 0 や 1 に近づくと、σ′(x) は 0 に近づく（つまり、勾配が消失する）ため、シグモイド関数は入力が大きく偏った場合に飽和しやすい特徴がある。

逆伝播において各層の誤差項は、前層の誤差項に対して、活性化関数の微分を掛け合わせる形で計算される。シグモイド関数を用いた場合、その微分は最大でも 0.25 であり、通常はそれ以下の値となる。従って、仮に各層での微分の平均値が 0.25 だとすると、層 L（最終層）からある層 l に伝播する際に、誤差項は理論上以下のように減衰する：

δ_l ∝ (0.25)^{(L − l)}

これは、各層で勾配に 0.25（またはそれ以下）のスケーリングがかかるため、層数が増えるとその積が指数関数的に小さくなり、特に入力層近くでは勾配がほとんど無視できるほど小さくなってしまうことを示している。

2. 現代における勾配消失問題の回避策

勾配消失問題は深層学習の大きな課題であった。しかし、様々な工夫により現在ではこの問題を軽減する手法が主流となっている。主な回避策は以下の通り：

• 活性化関数の選択：
  シグモイドやハイパボリックタンジェント（tanh）は出力が飽和しやすいのに対し、ReLU（Rectified Linear Unit）やその変種（Leaky ReLU, Parametric ReLU, ELU など）は、入力が正の場合ほぼ線形なため勾配が消失しにくくなっている。
  ReLU は以下のように定義される：
  

  ReLU(x) = max(0, x)

• 重みの初期化：
  Xavier初期化やHe初期化など、層毎の出力の分散を一定に保つ初期化手法を用いると、信号の前向き伝播や逆伝播中の値の極端な変動を抑え、勾配消失を防ぐ効果がある。
• 正規化技術：
  バッチ正規化（Batch Normalization）は、各ミニバッチごとに内部表現を正規化することで、各層の入力分布のばらつきを抑え、勾配消失や勾配爆発のリスクを低減する。
• 残差学習：
  ResNetではスキップ接続（shortcut connections）を導入することにより、勾配が直接初期層まで伝播できるため、非常に深いネットワークでも効果的に学習が進む。

誤差逆伝播法と自動微分

次は Back-Propagation です。同じくGPTさんに活躍してもらいました。計算もずいぶん正解するようになりました。

誤差逆伝播法は損失関数 $E$ の出力から各層のパラメータ（重みおよびバイアス）への偏微分（勾配）を効率よく求めるためのアルゴリズムで、DNNの基礎となっています。順伝播によりネットワークの出力が計算された後、損失 $E$に連鎖律（チェーンルール）を適用しながら逆方向に伝播させて各パラメータの更新に必要な勾配を求め、出力誤差を減少させる方向にパラメータを修正しながら学習を進めます。

 自動微分はフレームワークに実装されている勾配を求めるためのアルゴリズム。手計算では煩雑になりますが、フレームワークを利用することでシームレスかつ正確に扱えるようになります。

1. ネットワークの構成と順伝播

まず、シンプルな2層（隠れ層と出力層）の全結合ニューラルネットワークを考える。

• 隠れ層: 入力 $\mathbf{x}$ に対し、 $$\mathbf{z}^{(1)} = \mathbf{W}^{(1)}\, \mathbf{x} + \mathbf{b}^{(1)}, \quad \mathbf{a}^{(1)} = f\bigl(\mathbf{z}^{(1)}\bigr)$$ ここで、$f$ は隠れ層の活性化関数（例：シグモイドやReLU）。
• 出力層: 隠れ層の出力 $\mathbf{a}^{(1)}$ を受け、 $$\mathbf{z}^{(2)} = \mathbf{W}^{(2)}\, \mathbf{a}^{(1)} + \mathbf{b}^{(2)}, \quad \mathbf{a}^{(2)} = g\bigl(\mathbf{z}^{(2)}\bigr)$$ ここで、$g$ は出力層の活性化関数（例：シグモイド）。

損失関数 $E$（例：平均二乗誤差やクロスエントロピー誤差）によって、予測値 $\mathbf{a}^{(2)}$ と正解 $\mathbf{y}$ との誤差が評価される。

2. 誤差逆伝播法の基本概念

誤差逆伝播法では、連鎖律（チェーンルール）を用いて、損失関数 $E$ の内部変数に対する偏微分を段階的に求める。
すなわち、任意の中間変数 $a$ と $z$ について、以下の通り表せる。 $$\frac{\partial E}{\partial z} = \frac{\partial E}{\partial a} \cdot \frac{\partial a}{\partial z}$$

局所誤差を表すデルタ $\delta = \frac{\partial E}{\partial z}$ は、各層における勾配計算の基礎となる重要な指標である。これは、出力 $a$ の微小な変化 $\frac{\partial a}{\partial z}$が損失 $E$ にどのような影響を与えるか$\frac{\partial E}{\partial a}$ を明らかにしており、パラメータの更新方向および大きさを決定するのに欠かせない情報となっている。この考え方により、各層ごとのパラメータ更新に必要な勾配を効率的に計算できる点が、誤差逆伝播法の大きな利点である。

3. 出力層の誤差項δ

出力層では、出力 $\mathbf{a}^{(2)} = g(\mathbf{z}^{(2)})$ と正解との差が損失関数 $E$ で評価される。
たとえば、シグモイド関数の場合は、

$$g(z)=\frac{1}{1+e^{-z}}, \quad g'(z)=g(z)(1-g(z))$$

連鎖律を用いると、出力層でのデルタは次のように定義される：

$$\delta^{(2)} = \frac{\partial E}{\partial \mathbf{z}^{(2)}} = \frac{\partial E}{\partial \mathbf{a}^{(2)}} \odot g'\bigl(\mathbf{z}^{(2)}\bigr)$$

ここで$\odot$はアダマール積を意味しする。

この形式は、誤差の各成分が出力 $a$ に対する感度と、活性化関数 $g$ の変化率 $g'(z)$ とで表現されており、局所的な誤差伝播に必要な情報となる。

4. 出力層の重み $\mathbf{W}^{(2)}$ に関する偏微分

出力層においては、線形結合の定義より $$\mathbf{z}^{(2)} = \mathbf{W}^{(2)}\, \mathbf{a}^{(1)} + \mathbf{b}^{(2)}$$ となる。各出力ユニット $i$ の入力は、

$$z^{(2)}_i = \sum_j W^{(2)}_{ij}\, a^{(1)}_j + b^{(2)}_i$$

連鎖律により、各重み $W^{(2)}_{ij}$ に対する偏微分は、

$$\frac{\partial E}{\partial W^{(2)}_{ij}} = \frac{\partial E}{\partial z^{(2)}_i} \cdot \frac{\partial z^{(2)}_i}{\partial W^{(2)}_{ij}} = \delta^{(2)}_i \, a^{(1)}_j$$

すなわち、行列表現では

$$\frac{\partial E}{\partial \mathbf{W}^{(2)}} = \delta^{(2)} \Bigl(\mathbf{a}^{(1)}\Bigr)^\top$$

5. 隠れ層への逆伝播の導出

出力層で計算されたデルタ $\delta^{(2)}$ は、隠れ層に伝播される。ここでは、隠れ層の各ユニット $j$ におけるデルタ $\delta^{(1)}_j$とする。

まず、出力層の線形結合は次のように定義される： $$z^{(2)}_i = \sum_j W^{(2)}_{ij}\, a^{(1)}_j + b^{(2)}_i$$ そのため、隠れ層ユニット $j$ の出力 $a^{(1)}_j$ が出力層の $z^{(2)}_i$ に与える影響は $$\frac{\partial z^{(2)}_i}{\partial a^{(1)}_j} = W^{(2)}_{ij}\quad \text{となる。}$$

さらに、隠れ層では $\mathbf{a}^{(1)} = f(\mathbf{z}^{(1)})$ としており、活性化関数 $f$ の微分は $f'(z^{(1)}_j)$ となる。
連鎖律を適用すると、隠れ層ユニット $j$ における損失の微分は以下のように表される：

$$\frac{\partial E}{\partial z^{(1)}_j} = \sum_i \frac{\partial E}{\partial z^{(2)}_i}\,\frac{\partial z^{(2)}_i}{\partial a^{(1)}_j}\,\frac{\partial a^{(1)}_j}{\partial z^{(1)}_j} = \Bigl(\sum_i W^{(2)}_{ij}\,\delta^{(2)}_i\Bigr) \cdot f'\bigl(z^{(1)}_j\bigr)$$

ここで、隠れ層のデルタ $\delta^{(1)}_j$ を定義すると、

$$\delta^{(1)}_j = \frac{\partial E}{\partial z^{(1)}_j} = \Biggl(\sum_i W^{(2)}_{ij}\,\delta^{(2)}_i\Biggr) \cdot f'\bigl(z^{(1)}_j\bigr)$$

行列表現では、全隠れ層のデルタは次のようになる：

$$\delta^{(1)} = \Bigl(\mathbf{W}^{(2)}\Bigr)^\top\,\delta^{(2)} \odot f'\bigl(\mathbf{z}^{(1)}\bigr)$$

この導出から、デルタは「損失 $E$ の出力 $a$ に対する偏微分」と「活性化関数の微分 $\frac{\partial a}{\partial z}$」の積として定義できるため、各層の誤差伝播に必須の情報となる。
さらに、このデルタ情報を利用して、各パラメータの更新（モデルの改善）が効率的に行われる。

6. 計算例とモデル改善の流れ

以下に、1入力・1隠れ・1出力の場合の計算例を示す。ここではシグモイド関数を活性化関数とし、損失関数を平均二乗誤差 $$E=\frac{1}{2}\left(a^{(2)}-y\right)^2$$ とする。

[ステップ1] 順伝播

入力: $x=1$

重み・バイアス:
入力～隠れ層: $$W^{(1)}=0.5,\; b^{(1)}=0$$ 隠れ層～出力層: $$W^{(2)}=0.8,\; b^{(2)}=0$$

隠れ層:
線形結合:

$$z^{(1)} = 0.5 \times 1 + 0 = 0.5$$

シグモイド関数:

$$a^{(1)} = \sigma(0.5) \approx 0.6225$$

出力層:
線形結合:

$$z^{(2)} = 0.8 \times 0.6225 + 0 \approx 0.498$$

出力（シグモイド関数）:

$$a^{(2)} = \sigma(0.498) \approx 0.622$$

真の正解 $y=1$ として、損失は

$$E = 0.5 × (0.622 – 1)^2 ≈ 0.0714$$

[ステップ2] 出力層での逆伝播

損失の出力 $a^{(2)}$ に対する偏微分：

$$\frac{\partial E}{\partial a^{(2)}} = a^{(2)} - y \approx 0.622 - 1 = -0.378$$

シグモイドの微分：

$$g'(z^{(2)}) = a^{(2)}(1-a^{(2)}) \approx 0.622 \times 0.378 \approx 0.235$$

よって、出力層のデルタ：

$$\delta^{(2)} = (-0.378) \times 0.235 \approx -0.0888$$

これを用いて、出力層の重み勾配は：

$$\frac{\partial E}{\partial W^{(2)}} = \delta^{(2)} \times a^{(1)} \approx -0.0888 \times 0.6225 \approx -0.0553$$

[ステップ3] 隠れ層での逆伝播とモデル改善

次に、出力層から隠れ層へデルタを伝播させる。

$$\delta^{(1)} = (W^{(2)})^T \times \delta^{(2)} \odot f'(z^{(1)})$$

数値例では、$W^{(2)}=0.8$ のため、

$$(W^{(2)})^T \times \delta^{(2)} = 0.8 \times (-0.0888) \approx -0.07104$$

隠れ層のシグモイド微分:

$$f'(z^{(1)}) = a^{(1)}(1-a^{(1)}) \approx 0.6225 \times 0.3775 \approx 0.235$$

よって、隠れ層のデルタ：

$$\delta^{(1)} \approx -0.07104 \times 0.235 \approx -0.01669$$

隠れ層の重み勾配（入力 $x=1$ を用いる）：

$$\frac{\partial E}{\partial W^{(1)}} = \delta^{(1)} \times x \approx -0.01669 \times 1 = -0.01669$$

ここで得られた各層の勾配（出力層： -0.0553、隠れ層： -0.01669）から、パラメータ更新（例えば、確率的勾配降下法）を行う。
例えば学習率 $\eta$ を 0.1 とすると、重みの更新は以下のようになる：

$$W^{(2)}_{\text{new}} = W^{(2)} - \eta \times \frac{\partial E}{\partial W^{(2)}}, \quad W^{(1)}_{\text{new}} = W^{(1)} - \eta \times \frac{\partial E}{\partial W^{(1)}}$$

$$\begin{aligned} W^{(2)}_{\text{new}} &= 0.8 - 0.1 \times (-0.0553) \approx 0.80553,\\[8pt] W^{(1)}_{\text{new}} &= 0.5 - 0.1 \times (-0.01669) \approx 0.50167 \end{aligned}$$

この重みWの更新により、モデルは出力誤差を減少させる方向へ学習が進む。
すなわち、デルタが求めた各層での「局所的な誤差指標」を基にして、ネットワーク全体の性能が改善されることになる。

7. 自動微分

自動微分（Automatic Differentiation, AD）は、ニューラルネットワークやその他の複雑な計算グラフにおいて、各中間変数の勾配（微分値）を自動的に計算する技術である。これにより、手動での連鎖律（チェーンルール）の適用から解放され、正確かつ効率的な勾配計算が実現される。

自動微分の利点

• 正確性： 数値微分の近似誤差や記号微分の複雑な手作業を避け、連鎖律を忠実に適用して正確な勾配計算を行う。
• 効率性： 膨大なパラメータや複雑な計算グラフを持つモデルでも高速に勾配を計算できるため、大規模なニューラルネットワークの学習に最適。

自動微分の実装例：PyTorch

現代のディープラーニングライブラリ（例：PyTorch、TensorFlowなど）は自動微分機能を内蔵している。以下に、PyTorchを使った簡単な全結合ネットワークの例を示す。

import torch
 
# 計算グラフの構築
W1 = torch.tensor([[0.5]], requires_grad=True)
b1 = torch.tensor([0.0], requires_grad=True)
W2 = torch.tensor([[0.8]], requires_grad=True)
b2 = torch.tensor([0.0], requires_grad=True)
 
# 順伝播
x = torch.tensor([1.0])
z1 = W1 @ x + b1         # 隠れ層の線形結合
a1 = torch.sigmoid(z1)   # 活性化関数（シグモイド）
z2 = W2 @ a1 + b2        # 出力層の線形結合
a2 = torch.sigmoid(z2)   # 最終出力
 
# 逆伝播
y = torch.tensor([1.0])
loss = 0.5 * (a2 - y) ** 2
loss.backward()  # 逆伝播により、自動的に勾配計算が実行される
 

参考文献

• Rumelhart, Hinton, and Williams, "Learning representations by back-propagating errors" (1986)
• Goodfellow, Bengio, and Courville, "Deep Learning" (MIT Press, 2016)

シグモイド関数

あらためてシグモイド関数 by GPT。数式を打ち込む必要がなく、とても楽です。 

 　

1. ロジット関数と逆関数：ロジスティック関数（シグモイド関数）

2値分類（バイナリ分類）では、対象のデータを2つのクラス、通常は 0 と 1 に分類する。

ロジスティック回帰はこのような分類問題に用いられる手法で、線形予測子の出力を確率に変換するために、ロジスティック関数を利用する。また、対数オッズ（ロジット）として

$$\text{logit}(p) = \ln\left(\frac{p}{1-p}\right)$$

と定義され、これの逆関数としてロジスティック関数が用いられる。 $$z = \ln\left(\frac{p}{1-p}\right)$$ $$e^z = \frac{p}{1-p}$$ $$e^z (1-p) = p$$ $$e^z - e^z p = p$$ $$e^z = p + e^z p = p(1+e^z)$$ $$p=\frac{e^z}{1+e^z}$$ $$p=\frac{1}{1+e^{-z}}$$

2. シグモイド関数の微分

シグモイド関数を以下のように定義する。

$$S(x)=\frac{1}{1+e^{-x}}$$

シグモイド関数 $S(x)$ を微分するため、次の形に変形する。

$$S(x)=\left(1+e^{-x}\right)^{-1}$$

ここで、外部関数 $f(u)=u^{-1}$ と内部関数 $u(x)=1+e^{-x}$ と置くと、連鎖律により

$$\frac{dS}{dx}=\frac{df}{du}\cdot\frac{du}{dx}$$

外側の関数の微分

$$f(u)=u^{-1}\quad\Rightarrow\quad\frac{df}{du}=-u^{-2}$$

内側の関数の微分

$$u(x)=1+e^{-x}\quad\Rightarrow\quad\frac{du}{dx}=-e^{-x}$$

連鎖律の適用

$$\frac{dS}{dx}=-\left(1+e^{-x}\right)^{-2}\cdot\left(-e^{-x}\right) =\frac{e^{-x}}{\left(1+e^{-x}\right)^2}.$$

3. 微分式の形の簡略化

シグモイド関数の定義を用いると、その補数は $$1-S(x)=\frac{e^{-x}}{1+e^{-x}}$$ となり、積は $$S(x)\left(1-S(x)\right)=\frac{e^{-x}}{\left(1+e^{-x}\right)^2}$$ となる。

したがって、シグモイド関数の微分は

$$\frac{dS}{dx}=S(x)\left(1-S(x)\right)$$

と自身のみで表現可能。

回帰直線の傾きの導出と予測値の補正

PyCaret で上位５つの回帰モデルをチューニングし、結果を確認。
いずれも、予測値を縦軸、実測値を横軸にプロットすると傾きが１未満になります。理由があるはずだと調べてみましたが、ほとんど書かれていないようです。で、GPT君に聞きながら整理しました。

以下、実測値 $y$ と予測値 $\hat{y}$ の関係から回帰直線の傾き $\gamma$ の導出、その傾きが１以下となる理由、及び予測値が実測値の全変動を再現できない場合の補正方法について書き残しておきます。

※傾きの導出は「学習物理学入門」に書かれています。高校で習ったかな。

1. 傾きの公式の導出

ここでは、中心化したデータを用いて回帰直線の傾き $\gamma$ を求めるプロセスを示す。
まず、各データ点 $i$ に対して、実測値 $y_i$ と予測値 $\hat{y}_i$ の中心化された値を以下のように定義する：

$$\tilde{y}_i = y_i - \bar{y},\quad \tilde{\hat{y}}_i = \hat{y}_i - \bar{y}$$

このとき、中心化されたデータに対するモデルは

$$\tilde{\hat{y}}_i = \gamma\, \tilde{y}_i + \varepsilon_i$$

と表現され、ここで $\gamma$ は回帰直線の傾き、$\varepsilon_i$ は各データ点 $i$ における残差を示す。 残差の二乗和（Sum of Squared Errors, SSE）は次のように定義される：

$$\text{SSE}(\gamma) = \sum_{i=1}^{n} \left(\tilde{\hat{y}}_i - \gamma\, \tilde{y}_i\right)^2$$

1-1. SSE を $\gamma$ で最小化する過程

SSE を最小にする $\gamma$ を求めるため、まず SSE を $\gamma$ の関数として展開する。
各項は

$$f(\gamma) = \left[\tilde{\hat{y}}_i - \gamma\, \tilde{y}_i\right]^2$$

と表される。ここで、内部関数

$$u(\gamma) = \tilde{\hat{y}}_i - \gamma\, \tilde{y}_i$$

と置くと、$f(\gamma)$ は $f(\gamma) = \bigl[u(\gamma)\bigr]^2$ と表せる。

チェーンルールを用いて $f(\gamma)$ の $\gamma$ による微分を行う。まず、外側の関数 $g(u) = u^2$ の微分は

$$g'(u) = 2u$$

次に、内側の関数

$$u(\gamma) = \tilde{\hat{y}}_i - \gamma\, \tilde{y}_i$$

を $\gamma$ で微分する。ここで $\tilde{\hat{y}}_i$ は $\gamma$ に依存しない定数なので、その微分は 0 となり、

$$\frac{d}{d\gamma} u(\gamma) = -\tilde{y}_i$$

チェーンルールにより、全体の微分は

$$\frac{d}{d\gamma} f(\gamma) = g'\bigl(u(\gamma)\bigr) \cdot u'(\gamma)$$

これにより、具体的には

$$\frac{d}{d\gamma}\left[\tilde{\hat{y}}_i - \gamma\, \tilde{y}_i\right]^2 = 2\bigl(\tilde{\hat{y}}_i - \gamma\, \tilde{y}_i\bigr) \cdot (-\tilde{y}_i)$$

整理すると、

$$\frac{d}{d\gamma}\left[\tilde{\hat{y}}_i - \gamma\, \tilde{y}_i\right]^2 = -2\,\tilde{y}_i\,\bigl(\tilde{\hat{y}}_i - \gamma\, \tilde{y}_i\bigr)$$

これが、各項 $\left[\tilde{\hat{y}}_i - \gamma\, \tilde{y}_i\right]^2$ の $\gamma$ による微分の導出過程である。
次に、この導出を全データ点 $i = 1, 2, \ldots, n$ に対して合計し、SSE の全体の微分を求める。

$$\frac{d}{d\gamma}\text{SSE}(\gamma) = \sum_{i=1}^{n} -2\,\tilde{y}_i\,\bigl(\tilde{\hat{y}}_i - \gamma\, \tilde{y}_i\bigr)$$

SSE の最小値を与えるためには、この導関数を 0 とおく必要がある。すなわち、

$$-2 \sum_{i=1}^{n} \tilde{y}_i\,\bigl(\tilde{\hat{y}}_i - \gamma\, \tilde{y}_i\bigr) = 0$$

両辺を $-2$ で割ると、次の条件が得られる：

$$\sum_{i=1}^{n} \tilde{y}_i\,\bigl(\tilde{\hat{y}}_i - \gamma\, \tilde{y}_i\bigr) = 0$$

これを展開すると：

$$\sum_{i=1}^{n} \tilde{y}_i\, \tilde{\hat{y}}_i - \gamma \sum_{i=1}^{n} \tilde{y}_i^2 = 0$$

$\gamma$ について解くと、以下の式が得られる：

$$\gamma = \frac{\sum_{i=1}^{n} \tilde{y}_i\, \tilde{\hat{y}}_i}{\sum_{i=1}^{n} \tilde{y}_i^2}$$

分散と共分散の定義、 $\operatorname{Cov}(y,\hat{y}) = \frac{1}{n}\sum_{i=1}^{n} \tilde{y}_i\, \tilde{\hat{y}}_i$ および $\operatorname{Var}(y) = \frac{1}{n}\sum_{i=1}^{n} \tilde{y}_i^2$ を用いると、最終的に

$$\gamma = \frac{\operatorname{Cov}(y,\hat{y})}{\operatorname{Var}(y)}$$

となる。

1-2. 直交性による公式の変形

回帰分析では、予測値 $\hat{y}$ と残差 $e = y - \hat{y}$ が直交する性質、すなわち $\operatorname{Cov}(\hat{y}, e) = 0$ が成り立つ。
実測値は予測値と残差の和として表せるので、

$$y = \hat{y} + e$$

となる。この関係から、共分散の線形性を用いると、

$$\operatorname{Cov}(y,\hat{y}) = \operatorname{Cov}(\hat{y}+e,\hat{y}) = \operatorname{Cov}(\hat{y},\hat{y}) + \operatorname{Cov}(e,\hat{y})$$

ここで $\operatorname{Cov}(e,\hat{y}) = 0$ であるため、

$$\operatorname{Cov}(y,\hat{y}) = \operatorname{Var}(\hat{y})$$

従って、もともとの $\gamma$ の式は

$$\gamma = \frac{\operatorname{Var}(\hat{y})}{\operatorname{Var}(y)}$$

と表される。

2. 予測値の補正方法とその導出

通常、モデルの予測値 $\hat{y}$ は実測値 $y$ の全てのばらつきを再現できない。その理由は、回帰モデルが実際のデータのばらつきを完全には捉えられず、残差 $e$ すなわちモデルでは説明できない成分が存在するためである。

実測値は予測値と残差の和として表現できる：

$$y = \hat{y} + e$$

ここで、回帰分析の基本的性質として、予測値 $\hat{y}$ と残差 $e$ は直交すると仮定され、$\operatorname{Cov}(\hat{y},e) = 0$ が成立する。
この性質を用いて、実測値 $y$ の分散を計算すると、

$$\operatorname{Var}(y) = \operatorname{Var}(\hat{y}+e) = \operatorname{Var}(\hat{y}) + \operatorname{Var}(e) + 2\,\operatorname{Cov}(\hat{y},e)$$

ここで $\operatorname{Cov}(\hat{y},e)=0$ なので、

$$\operatorname{Var}(y) = \operatorname{Var}(\hat{y}) + \operatorname{Var}(e)$$

残差の分散 $\operatorname{Var}(e)$ は非負のため、

$$\operatorname{Var}(\hat{y}) \le \operatorname{Var}(y)$$

すなわち、モデルの予測値 $\hat{y}$ は実測値 $y$ の全ばらつきを再現できず、通常は $\operatorname{Var}(\hat{y}) < \operatorname{Var}(y)$ となる。
この状態では、中心化された回帰直線の傾きは

$$\gamma = \frac{\operatorname{Var}(\hat{y})}{\operatorname{Var}(y)} < 1$$

この差を補正するために、予測値の分散を実測値の分散に合わせる補正が必要となる。

2-1. 標準偏差比による補正

実測値 $y$ の標準偏差を $\sigma_y$、予測値 $\hat{y}$ の標準偏差を $\sigma_{\hat{y}}$ とすると、
補正後の予測値は以下の式で定義される：

$$\hat{y}_{\text{new}} = \bar{y} + \Bigl(\hat{y} - \bar{y}\Bigr) \times \frac{\sigma_y}{\sigma_{\hat{y}}}$$

2-2. 補正方法の導出（標準偏差比による補正の意味）

まず、予測値 $\hat{y}$ を中心化すると $\hat{y} - \bar{y}$ となる。
補正前の分散は $\operatorname{Var}(\hat{y}) = \sigma_{\hat{y}}^2$ である。
補正係数 $\frac{\sigma_y}{\sigma_{\hat{y}}}$ を用いることで、補正後の分散は

$$\operatorname{Var}(\hat{y}_{\text{new}}) = \left(\frac{\sigma_y}{\sigma_{\hat{y}}}\right)^2 \operatorname{Var}(\hat{y}) = \left(\frac{\sigma_y}{\sigma_{\hat{y}}}\right)^2 \sigma_{\hat{y}}^2 = \sigma_y^2$$

これにより、補正後の予測値は実測値 $y$ の全体のばらつきに一致し、中心化された回帰直線の傾きが 1 に調整される。

3. 参考情報

3-1. 相関係数と傾きの違い

2つの変数 $y$ と $\hat{y}$ の相関係数は通常、

$$r_{y,\hat{y}} = \frac{\operatorname{Cov}(y,\hat{y})}{\sigma_y\, \sigma_{\hat{y}}}$$

と定義される。一方、傾きは実測値 $y$ の分散を基準としているため、

$$\gamma = \frac{\operatorname{Cov}(y,\hat{y})}{\sigma_y^2}$$

と表現され、尺度が異なるため相関係数とは意味が異なる。

3-2. 決定係数とのつながり

決定係数 $R^2$ は、モデルが実測値のばらつきをどの程度説明できるかを示す指標であり、以下のように定義される：

$$R^2 = \frac{\operatorname{Var}(\hat{y})}{\operatorname{Var}(y)}$$

未補正の場合は、回帰直線の傾き $\gamma$ と同様の形となり、通常 $R^2 < 1$ となる。

4. まとめ

・中心化したデータに対して最小二乗法を用い、残差の二乗和（SSE）を最小にすることで回帰直線の傾きは

$$\gamma = \frac{\operatorname{Cov}(y,\hat{y})}{\operatorname{Var}(y)}$$

と導出される。さらに、残差と予測値の直交性によりこの式は

$$\gamma = \frac{\operatorname{Var}(\hat{y})}{\operatorname{Var}(y)}$$

とも表現される。
通常、モデルは実測値の全変動を再現できないため、$\operatorname{Var}(\hat{y}) < \operatorname{Var}(y)$ となり、補正前の $\gamma$ は 1 未満となる。

・補正方法の一つとして、標準偏差比を用いて予測値の分散を実測値の分散に合わせる方法があり、

$$\hat{y}_{\text{new}} = \bar{y} + \Bigl(\hat{y} - \bar{y}\Bigr) \times \frac{\sigma_y}{\sigma_{\hat{y}}}$$

により、補正後の予測値は実測値 $y$ のばらつき（$\sigma_y^2$）と一致し、中心化された回帰直線の傾きが理想的な 1 に調整される。

ラズパイの時刻設定

Raspberry Pi 3B/4 にて定期的な時計合わせを行いたいと、GPT 君に教えを請いました。

• 現在はインターネットにつながっており、chronyを使用中
• ネットにつながらない時は GPS モジュールと gpsd + chrony を使う予定
• RTCを組み込んで、起動時にRTCから本体の時刻を修正、定期的にRTCの時刻を修正
• GPIOは他のモジュールに占用されており、RTC、GPSはUSB接続

このような条件です。

※RTC：東京デバイセズ TDPC0205

※WinSCP
ログイン画面の[設定]
[環境]->[シェル]
「シェル」が[デフォルト]から[sudo su -]に変更

1. GPS

1-1. update_time.sh のコピー

• script フォルダをデスクトップに作成、書き込みパーミッションを与える
• update_time.sh（及び test.sh）をコピー、実行パーミッションを与える
• update_time.sh には、起動時に RTC からシステム時刻に書き込む設定、11分毎にシステム時刻をRTCへ書き込む設定(Chronyのrtcsyncに対応していないRTCの場合)を記載

1-2. gpsd および gpsd-clients の導入

• USB 接続の GPS モジュール使用のため、以下の関連パッケージをインストール

sudo apt install gpsd gpsd-clients

1-3. gpsd の設定

• VID/PIDの調査

lsusb→”aaaa:bbbb”

• /etc/udev/rules.d/100-usb-GPS_RTC.rules ファイルを作成し、次の内容を記載

SUBSYSTEM=="tty", ATTRS{idVendor}=="aaaa", ATTRS{idProduct}=="bbbb", SYMLINK+="gps0"

• 再起動後、NMEA 文（$GPGGA, $GPRMC 等）が表示されるか確認

sudo cat /dev/gps0

• GPS モジュールの接続先に合わせ、設定ファイルを編集
  /etc/default/gpsd の例

START_DAEMON="true"
GPSD_OPTIONS="-b -n"
DEVICES="/dev/gps0"
USBAUTO="false"

• 設定後、サービスを再起動。

sudo systemctl restart gpsd

1-4. SHMから読み取る設定を追加（gpsdから取得）

• /etc/chrony/chrony.confに追記

# SHM(shared memory)から読み取る設定
refclock SHM 0 refid GPS

• 設定後、サービスの再起動を行う

sudo systemctl restart chrony chronyd

2.RTC

2-1. 東京デバイセズ制御プログラムのインストール

• https://github.com/tokyodevices/td-usb/blob/master/README_ja.md

sudo apt install libusb-dev

• scriptフォルダからターミナルを起動

git clone https://github.com/tokyodevices/td-usb
cd td-usb
make

2-2. 100-usb-GPS_RTC.rulesに追記

• VID/PIDの調査

dmesg | grep TDPC0205

もしくは lsusb→”cccc:dddd”

• /etc/udev/rules.d/100-usb-GPS_RTC.rules ファイル

SUBSYSTEM=="usb", ATTRS{idVendor}=="cccc", ATTRS{idProduct}=="dddd", MODE="0666"

2-3. 起動時：RTCから時刻設定

• /etc/systemd/system/time_sync_RTC.serviceファイルを作成し、次の内容を記載

[Unit]
Description=起動時にRTCの時刻をシステム時刻に反映するサービス
After=network.target

[Service]
Type=oneshot
ExecStartPre=/bin/sleep 20
ExecStart=/home/pi/Desktop/script/update_time.sh boot

[Install]
WantedBy=multi-user.target

2-4. 有効化＆起動

sudo systemctl daemon-reload
sudo systemctl enable time_sync_RTC.service

3. systemd タイマーによる定期実行

• 11 分毎に時刻同期を行う場合、systemd の service と timer を利用して自動実行する

3-1. /etc/systemd/system/time_sync.service を作成

[Unit]
Description=時刻同期サービス

[Service]
Type=oneshot
ExecStart=/home/pi/Desktop/script/update_time.sh

3-2. /etc/systemd/system/time_sync.timer を作成

[Unit]
Description=11分ごとに時刻同期を実施するタイマー

[Timer]
OnBootSec=1min
OnUnitActiveSec=11min

[Install]
WantedBy=timers.target

3-3. タイマーを有効化＆起動

sudo systemctl daemon-reload
sudo systemctl restart time_sync.service
sudo systemctl status time_sync.service

• time_sync.timerの起動状況確認

sudo systemctl status time_sync.timer

※時刻関連の確認

• Chrony は ntpd のような「fudge」コマンドによる時刻オフセット補正をサポートしていない
• GPS が UTC で提供する時刻とシステム時刻（内部は通常 UTC）の整合性をシステムタイムゾーン設定で管理
• システムのタイムゾーン確認

timedatectl

• タイムゾーンを設定

sudo timedatectl set-timezone Asia/Tokyo

• 設定変更後は、Chrony デーモンを再起動して反映

sudo systemctl restart chronyd

• 接続確認

chronyc sources -v

***********************************************
20250326 RTC部分を追加、全体を修正

データ解釈の基礎

データを把握するために、統計的手法をよく利用しています。

もともと統計が好きでなく（どちらかというと避けてきたので）、そちらの基礎知識が不足していました。そのため、図書を読みながら試行を重ね、手順を探りながら経験を積んできました。データ分析や機械学習に関しても同様です。基本的には独習が中心となるため、どうしても抜けのある知識になります。

昨日、本屋で以下の2冊を購入しました。

指標・特徴量の設計から始めるデータ可視化学入門　データを洞察につなげる技術（江崎貴裕） | ソシム
データ分析のための数理モデル入門本質をとらえた分析のために（東京大学先端科学技術研究センター 江崎貴裕）| ソシム

入門書ですので"広く浅く"といった内容ですが、抜けをチェックするのに都合の良い図書でした。見落としが少なく安心した一方で、こんな簡単な点や用語をこぼしていた！ということも。
いくつか補間できましたので、次に進みましょう。

崩壊土砂量の推定式

他支店の方から推定崩壊土砂量データを頂きました。

分布を図化すると、少なからず異常値が含まれていそうでした。式の適用範囲を聞いても知らないとのことでしたので、こちらで調べることにしました。

土砂量式は時々見かける以下の形です。αとγはフィッティングパラメータです。原著は日本人で、1969年発表ですから50年以上前の古い式です。豪雨時の崩壊がソースで、α=0.83, γ=1.18 です。面積は 500～100,000m2 の範囲です。
Frequency Distribution of the Magnitude of the Landslides Caused by Heavy Rain-Fall
$$\begin{align} V=\alpha\ A^\gamma \end{align}$$

他国でもよく利用されているようです。いくつかの文献で求められた値が modified from Guzzetti et al. 2009 として以下に掲載されています。
Full article: Landslides detection and volume estimation in Jinbu area of Korea

この中の Guzzetti et al. 2009 の値、α=0.074, γ=1.45が例として 深層崩壊の発生する恐れのある斜面抽出技術手法及びリスク評価手法に関する研究－土木研究所 に掲載されています。おそらく、ソースの範囲 2～１,000,000,000m2 が広いからでしょう。 

「たとえば」として書かれていても、土研資料に取り上げられた文献値に頼りたくなるのでしょう。今回もこの式が利用されていました。範囲としては問題なし。うーん。

H/V の深度表現 その2

先月計算していたH/Vの深度情報を可視化してみました。以下の簡易手法の結果です。
https://phreeqc.blogspot.com/2025/02/hv.html

1つ目。
流れ版のような構造です。面の真ん中と奥が少しへこんで複数の段差が見えます。想定していた地すべりブロックの底面に対応しています。が、少し浅い。すべった後に動きがなかった現場ですので深めに想定していたのですが、実際はこの程度の深度なのかもしれません。
この結果を使ってすべり面を設定するのは少し難しいと思います。

２つ目
こちらも１つの面の中で段差が生じています。

この現場は２方向への動きがあったようですが、それに対応しているようです。
すべり面の中～下部に対応する面が現れているようですが、頭部は見えません。傾斜がきついとH/Vのピークがうまく出てこないのでしょう。

参考

・10度以上で不明瞭

・基盤層と堆積層のS波速度コントラストが概ね3倍程度以上。

https://phreeqc.blogspot.com/2017/04/hv.html

地表移動ベクトルを用いた3次元すべり面補正

地すべり縦断に対し、SLBLですべりラインを作成したのちに地表の観測結果で勾配を修正するコードを見かけました。

後半の地表変位を利用してすべり面を作成する方法は、昔、土研さんが研究成果を出版されています。プログラムも公開されています。残念ながら、地すべり事業で地表面変位を測ることが標準とされていなかったため、当時はお金をかけてデータをとることが一般的でなく、活用する場に恵まれませんでした。
すべり線推定システム｜地すべり 土木研究所

近年、まだ制限はあるものの３次元で広域に地表面変位量を求めることが可能となりました。SAR や LiDAR ＋ PIV で移動ベクトルを算出することが可能です。これらを利用してすべり面を推定している文献を国外では見かけます。

SLBL や GEORAMA ですべり面サーフェスを作成した後に、勾配等で微修正するだけなら容易でしょう。
ということで、ペナルティ付き4項を足して最小化するだけの単純なスクリプトをGPTさんに作ってもらいました。

正則化項    Reg ：標高変化は小さいほうが良い
滑らかさ項 S：急激な変化が小さいほど良い
標高項       H ：指定標高に近いほど良い
勾配項       G：指定勾配に近いほど良い

滑らかさの指標をどうするかで迷いましたが、①4方向ラプラシアンと②ガウス関数でターゲット周辺に重みを付けて修正した場合の2パターンを試しました。

オリジナル

①ラプラシアン

②重み付き評価

どちらの補正方法が良いのかわかりません。が、意図したように動いています。

実際は3次元移動場がすべり面傾斜方向に制約を受けているのかどうかは不明です。数を試すしかないでしょうね。すべり面と移動量を扱うことがあれば試してみましょう。

Landslide ＋ 機械学習の傾向

Landslide identification using machine learning techniques: Review, motivation, and future prospects | Earth Science Informatics

大きく分けて、２つの目的で機械学習が使われていました。

• 発生場所の特定（Landslide detection, identification, mapping）
• 危険度予測（Landslide susceptibility prediction, mapping）

2022年までですので、画像を扱う場合はCNNが多いようです。機械学習分野で利用されているの ViT ベースのアーキテクチャは見当たりません。

Landslide susceptibility は2003年から。2017年までは ANN と SVM が主体。その後はCNN, RF, DNN, GBDTなど。現在は LightGBM でしょうか。 

少し古い情報ですが、変遷の一例として使えそうです。

Year	Objective	Machine Learning Model
2003	Landslide susceptibility analysis	GIS and ANN
2005	Landslide susceptibility analysis	Multiple ANN and GIS
2005	Landslide susceptibility assessment	ANN
2005	Landslide susceptibility mapping	Neural networks and GIS
2006	Landslide susceptibility analysis	SVM
2007	Landslide susceptibility analysis	ANN and GIS
2008	Landslide susceptibility analysis	SVM
2008	Landslide susceptibility analysis and prediction	GMM
2008	Landslide susceptibility zonation	ANN and DA
2009	Landslide susceptibility analysis	Decision tree
2011	Landslide susceptibility zonation	SVM
2012	Landslide susceptibility assessment	SVM, Decision Tree, and Logistic
2013	Landslide susceptibility mapping	SVM
2014	Landslide susceptibility mapping	Decision tree
2014	Landslide susceptibility mapping	Decision tree
2016	Landslide susceptibility mapping	ANN
2017	Landslide susceptibility assessment	SVM
2017	Landslide susceptibility mapping	SVM
2018	Landslide susceptibility mapping	Genetic Algorithm
2019	Landslide susceptibility mapping	Random forest, Naive Bayes, and Boosted Logistic Regression
2019	Landslide susceptibility mapping	CNN
2020	Landslide susceptibility prediction	FC-SAE
2021	Landslide susceptibility mapping	Random forest
2021	Landslide susceptibility prediction	GB
2021	Landslide susceptibility prediction	CNN
2022	Landslide susceptibility mapping	CNN and GWO
2022	Landslide susceptibility mapping	Random Forest
2022	Landslide susceptibility mapping	Deep Neural Networks
2018a	Landslide susceptibility assessment	Random forest
2018b	Landslide susceptibility assessment	Boosted tree