深層学習で生成したデータはベンフォードの法則を突破できるのか

みなさんが大好きなベンフォードの法則を使用して深層学習で生成したデータがこれに則るのか検証してみたいと思います。

## ベンフォードの法則
ベンフォードの法則は、自然界における数値データにおいて、先頭の桁（1から9）が現れる頻度についての統計的な法則です。数値データの先頭の桁は1である確率が最も高く、9である確率が最も低くなります。

例えば、人口統計、株価、物理定数など、様々なデータがこの法則に当てはまります。

## 深層学習による画像生成
深層学習を用いることで、実在しない画像を生成することが可能になりました。画像生成には、GAN (敵対的生成ネットワーク)と呼ばれる技術がよく使われます。 GANは、GeneratorとDiscriminatorという2つのネットワークから構成されます。
Generator：ランダムなノイズから画像を生成
Discriminator：入力された画像が本物か偽物かを判定
GeneratorとDiscriminatorは互いに競い合うように学習するため、最終的にGeneratorは本物と見分けがつかないような高精度な画像を生成できるようになります。

## 実験：深層学習モデル vs ベンフォードの法則
今回は、深層学習モデルにDCGANを用い、CIFAR-10データセットで学習を行います。そして、生成された画像データと、オリジナルのCIFAR-10データセットそれぞれに対して、ベンフォードの法則適合度p値を評価します。

p値は、帰無仮説が正しいという前提のもとで、観測されたデータよりも極端なデータが得られる確率を表し、今回は有意水準を慣習的な0.05として統計的に判断します。

プログラムは以下通りです。

```python
import torch
import torch.nn as nn
import torch.optim as optim
from torchvision import datasets, transforms
from torch.utils.data import DataLoader
from collections import Counter
import math
from scipy.stats import chisquare

# ベンフォードの法則の理論値
benford_distribution = [math.log10(1 + 1 / d) for d in range(1, 10)]

# データの最上位桁の出現頻度を計算する関数
def calculate_first_digit_frequency(data):
  first_digits = [int(str(abs(x))[0]) for x in data if x != 0]
  return list(Counter(first_digits).values())

# ベンフォードの法則適合度を評価する関数
def evaluate_benford_fit(data):
  observed_frequencies = calculate_first_digit_frequency(data)
  expected_frequencies = [len(data) * p for p in benford_distribution]
  _, p_value = chisquare(observed_frequencies, f_exp=expected_frequencies)
  return p_value

# DCGANの定義
class Generator(nn.Module):
    def __init__(self, nz, ngf, nc):
        super(Generator, self).__init__()
        self.main = nn.Sequential(
            # 入力はZで畳み込みに入る
            nn.ConvTranspose2d( nz, ngf * 8, 4, 1, 0, bias=False),
            nn.BatchNorm2d(ngf * 8),
            nn.ReLU(True),
            # 状態サイズ: (ngf*8) x 4 x 4
            nn.ConvTranspose2d(ngf * 8, ngf * 4, 4, 2, 1, bias=False),
            nn.BatchNorm2d(ngf * 4),
            nn.ReLU(True),
            # 状態サイズ: (ngf*4) x 8 x 8
            nn.ConvTranspose2d( ngf * 4, ngf * 2, 4, 2, 1, bias=False),
            nn.BatchNorm2d(ngf * 2),
            nn.ReLU(True),
            # 状態サイズ: (ngf*2) x 16 x 16
            nn.ConvTranspose2d( ngf * 2, ngf, 4, 2, 1, bias=False),
            nn.BatchNorm2d(ngf),
            nn.ReLU(True),
            # 状態サイズ: (ngf) x 32 x 32
            nn.ConvTranspose2d( ngf, nc, 4, 2, 1, bias=False),
            nn.Tanh()
            # 状態サイズ: (nc) x 64 x 64
        )

    def forward(self, input):
        return self.main(input)

class Discriminator(nn.Module):
    def __init__(self, nc, ndf):
        super(Discriminator, self).__init__()
        self.main = nn.Sequential(
            # 入力は (nc) x 64 x 64
            nn.Conv2d(nc, ndf, 4, 2, 1, bias=False),
            nn.LeakyReLU(0.2, inplace=True),
            # 状態サイズ: (ndf) x 32 x 32
            nn.Conv2d(ndf, ndf * 2, 4, 2, 1, bias=False),
            nn.BatchNorm2d(ndf * 2),
            nn.LeakyReLU(0.2, inplace=True),
            # 状態サイズ: (ndf2) x 16 x 16
            nn.Conv2d(ndf * 2, ndf * 4, 4, 2, 1, bias=False),
            nn.BatchNorm2d(ndf * 4),
            nn.LeakyReLU(0.2, inplace=True),
            # 状態サイズ: (ndf4) x 8 x 8
            nn.Conv2d(ndf * 4, ndf * 8, 4, 2, 1, bias=False),
            nn.BatchNorm2d(ndf * 8),
            nn.LeakyReLU(0.2, inplace=True),
            # 状態サイズ: (ndf*8) x 4 x 4
            nn.Conv2d(ndf * 8, 1, 4, 1, 0, bias=False),
            nn.Sigmoid()
        )

    def forward(self, input):
        return self.main(input).view(-1, 1).squeeze(1)

# ハイパーパラメータの設定
batch_size = 64
image_size = 64
nc = 3 # カラー画像なのでチャネル数は3
nz = 100 # 入力ノイズの次元数
ngf = 64 # Generatorのフィルター数
ndf = 64 # Discriminatorのフィルター数
num_epochs = 100 # 学習回数
lr = 0.0002 # 学習率
beta1 = 0.5 # Adamオプティマイザにおける指数移動平均を直近の勾配に依存するかの度合い

# デバイスの設定
device = torch.device("cuda" if torch.cuda.is_available() else "cpu")

# データセットの読み込みと前処理
transform = transforms.Compose([
    transforms.Resize(image_size),
    transforms.ToTensor(),
    transforms.Normalize((0.5, 0.5, 0.5), (0.5, 0.5, 0.5))
])

train_dataset = datasets.CIFAR10(root='./data', train=True, download=True, transform=transform)
train_loader = DataLoader(train_dataset, batch_size=batch_size, shuffle=True)

# モデルのインスタンス化
netG = Generator(nz, ngf, nc).to(device)
netD = Discriminator(nc, ndf).to(device)

# 最適化関数の定義
optimizerD = optim.Adam(netD.parameters(), lr=lr, betas=(beta1, 0.999))
optimizerG = optim.Adam(netG.parameters(), lr=lr, betas=(beta1, 0.999))

# 損失関数の定義
criterion = nn.BCELoss()

# 学習
for epoch in range(num_epochs):
    for i, (real_images, _) in enumerate(train_loader):
        # データをデバイスに転送
        real_images = real_images.to(device)
        batch_size = real_images.size(0)

        # Discriminatorの学習
        netD.zero_grad()
        label = torch.full((batch_size,), 1, dtype=torch.float, device=device)
        output = netD(real_images)
        errD_real = criterion(output, label)
        errD_real.backward()

        noise = torch.randn(batch_size, nz, 1, 1, device=device)
        fake_images = netG(noise)
        label.fill_(0)
        output = netD(fake_images.detach())
        errD_fake = criterion(output, label)
        errD_fake.backward()
        errD = errD_real + errD_fake
        optimizerD.step()

        # Generatorの学習
        netG.zero_grad()
        label.fill_(1)
        output = netD(fake_images)
        errG = criterion(output, label)
        errG.backward()
        optimizerG.step()

        # 進捗表示
        if i % 100 == 0:
            print('[%d/%d][%d/%d] Loss_D: %.4f Loss_G: %.4f'
                  % (epoch + 1, num_epochs, i, len(train_loader),
                     errD.item(), errG.item()))

# 生成された画像データの評価
generated_images = netG(torch.randn(10000, nz, 1, 1, device=device)).detach().cpu()
generated_data = generated_images.flatten().tolist() 
p_value_generated = evaluate_benford_fit(generated_data)

# 自然画像データセットの評価
natural_data = train_dataset.data.flatten().tolist()
p_value_natural = evaluate_benford_fit(natural_data)

print(f'自然画像データセットのベンフォードの法則 p値: {p_value_natural}')
print(f'生成された画像データのベンフォードの法則 p値: {p_value_generated}')
```

## 結果
自然画像のp値は0.02562
生成画像のp値は0.15124
自然画像は有意水準内であるが、生成画像は有意水準を超える結果となりました。

## 考察
生成画像のp値が有意水準の範囲を超えた理由として、生成モデルがまだ十分に学習できていない、または、複雑な分布を生成する能力が不足していることが考えられます。また、使用している損失関数が、ベンフォードの法則のような特定の分布への適合性を直接的に最適化していない可能性も考えられます。
加えて、深層学習モデルは学習データの分布を再現するように学習するため、学習データに存在しないパターンを生成することは苦手です。ベンフォードの法則は、自然界のデータに広く見られるパターンですが、深層学習モデルはこの法則を学習していないため、生成データに反映されない可能性があるため、このような結果になったと考えられます。

## まとめ
今回の評価では、深層学習で生成されたデータが、必ずしも自然なデータと同じ性質を持つわけではないことが示唆されました。ただし、これはあくまでも統計的検定になるので、別の評価方法や有意水準の設定を試してみると違った結果が出て面白いかもしれません。