乙Py先生のプログラミング教室
初学者のためのプログラミング学習サイト

プログラミング

今回もJSでクソゲーを作ってみました


解説動画はこちら







クソゲー作ってみた

その名も
「ずっこんテトリス」
です

普通のテトリスを改変して
少し面白くなるようにしてみました
こちらからプレイできます。

ずっこんテトリス

上下左右キーの操作のみの
シンプル設計で

長いテトリス棒を追加しました
sample


こやつは「ずっこんテトリス棒」ですね

コードはこんな感じになっています。
const game = document.getElementById("game");
const rareChanceInput = document.getElementById("rareChance");
const chanceValueDisplay = document.getElementById("chanceValue");
const ROWS = 20;
const COLS = 10;

// グリッドを初期化
const grid = Array.from({ length: ROWS }, () => Array(COLS).fill(0));

// 通常のテトリミノ
const tetrominoes = [
  [[1, 1, 1, 1]], // I
  [
    [1, 1],
    [1, 1],
  ], // O
  [
    [0, 1, 0],
    [1, 1, 1],
  ], // T
  [
    [1, 1, 0],
    [0, 1, 1],
  ], // S
  [
    [0, 1, 1],
    [1, 1, 0],
  ], // Z
  [
    [1, 1, 1],
    [1, 0, 0],
  ], // L
  [
    [1, 1, 1],
    [0, 0, 1],
  ], // J
];

// レアな「ロング棒」長さを8に変更
const rareTetromino = [
  [[1, 1, 1, 1, 1, 1, 1, 1]], // 長さ8の棒
];

// テトリミノの状態
let currentTetromino = getRandomTetromino();
let currentRow = 0;
let currentCol = Math.floor((COLS - currentTetromino[0].length) / 2);

// グリッドを描画
function drawGrid() {
  game.innerHTML = "";
  for (let row = 0; row < ROWS; row++) {
    for (let col = 0; col < COLS; col++) {
      const cell = document.createElement("div");
      cell.classList.add("cell");
      if (grid[row][col] === 1) {
        cell.classList.add("filled");
      }
      game.appendChild(cell);
    }
  }
}

// テトリミノを描画
function drawTetromino() {
  currentTetromino.forEach((row, r) => {
    row.forEach((value, c) => {
      if (value && currentRow + r >= 0) {
        grid[currentRow + r][currentCol + c] = 1;
      }
    });
  });
}

// テトリミノを削除
function clearTetromino() {
  currentTetromino.forEach((row, r) => {
    row.forEach((value, c) => {
      if (value && currentRow + r >= 0) {
        grid[currentRow + r][currentCol + c] = 0;
      }
    });
  });
}

// 衝突判定
function isValidMove(newRow, newCol, newTetromino) {
  return newTetromino.every((row, r) =>
    row.every((value, c) => {
      const x = newCol + c;
      const y = newRow + r;
      return (
        !value ||
        (y >= 0 && y < ROWS && x >= 0 && x < COLS && grid[y][x] === 0)
      );
    })
  );
}

// ラインを削除
function clearLines() {
  for (let row = ROWS - 1; row >= 0; row--) {
    if (grid[row].every((cell) => cell === 1)) {
      grid.splice(row, 1);
      grid.unshift(Array(COLS).fill(0));
      row++;
    }
  }
}

// ランダムなテトリミノを取得(レア形状の低確率出現を含む)
function getRandomTetromino() {
  const rareChance = parseFloat(rareChanceInput.value); // スライダーの値を取得
  if (Math.random() < rareChance) { // ロング棒の出現確率
    return rareTetromino[0];
  } else {
    return tetrominoes[Math.floor(Math.random() * tetrominoes.length)];
  }
}

// テトリミノを回転
function rotateTetromino() {
  const newTetromino = currentTetromino[0].map((_, colIndex) =>
    currentTetromino.map((row) => row[colIndex]).reverse()
  );

  if (isValidMove(currentRow, currentCol, newTetromino)) {
    currentTetromino = newTetromino;
  }
}

// ゲームのループ
function gameLoop() {
  clearTetromino();
  if (isValidMove(currentRow + 1, currentCol, currentTetromino)) {
    currentRow++;
  } else {
    drawTetromino();
    clearLines();

    // 次のテトリミノを生成
    currentTetromino = getRandomTetromino();
    currentRow = 0;
    currentCol = Math.floor((COLS - currentTetromino[0].length) / 2);

    // ゲームオーバー判定
    if (!isValidMove(currentRow, currentCol, currentTetromino)) {
      alert("Game Over");
      grid.forEach((row) => row.fill(0));
      currentTetromino = getRandomTetromino();
      currentRow = 0;
      currentCol = Math.floor((COLS - currentTetromino[0].length) / 2);
    }
  }
  drawTetromino();
  drawGrid();
}

// キー操作
document.addEventListener("keydown", (e) => {
  clearTetromino();
  if (e.key === "ArrowLeft" && isValidMove(currentRow, currentCol - 1, currentTetromino)) {
    currentCol--;
  } else if (e.key === "ArrowRight" && isValidMove(currentRow, currentCol + 1, currentTetromino)) {
    currentCol++;
  } else if (e.key === "ArrowDown") {
    if (isValidMove(currentRow + 1, currentCol, currentTetromino)) {
      currentRow++;
    }
  } else if (e.key === "ArrowUp") {
    rotateTetromino();
  }
  drawTetromino();
  drawGrid();
});

// スライダーの値を表示
rareChanceInput.addEventListener("input", () => {
  chanceValueDisplay.textContent = parseFloat(rareChanceInput.value).toFixed(2);
});

// ゲーム開始
setInterval(gameLoop, 500);
drawGrid();


ここでレアなテトリス棒として定義し
出現確率を画面上のスライダーの値で設定しています。

ここを上げると、出現確率も上がります。

どの確率が一番面白くなるかは分かりませんので
色々変えて試していただければと思います。

コードも改変すると
面白くなるかもしれないので
色々遊んでみてください

それでは






 

今回は簡単なゲームを
Javascriptで作ってみました

解説動画はこちら




今回は「ナンスピ」と呼ばれるゲームを再現してみました。

4972825233115-3-1400x1400のコピー


東急ハンズとかで売られている
トイグッズになりますが

1から25までの数字をぽちぽち押していって
タイムを競うものです

製品のHPはこちら

ナンスピ

試し書きには
20秒以下 : げべげべ
20秒以内 : 平均的
16秒以内 : 良い
12秒以内 : すごく良い
6.25秒 : パイロットの平均

とありました。

買うと高いので
Javascriptのコードにしてみました


コードはこちら
(Codepen)

ソースコード

HTML,CSS
Javascriptを一緒くたにしています。

ソースコードをコピーして
HTMLファイル保存してブラウザ実行するか
Codepenでそのまま楽しめます。

スクリーンショット 2025-01-18 16.15.58


単純に「ゲーム開始」ボタンで
ゲームスタート

1から25までの数字をクリックしてください

スクリーンショット 2025-01-18 16.16.09


25までクリックし終えるとタイムがポップアップします。
スクリーンショット 2025-01-18 16.15.34
過去の記録もされるので
過去の自分を超えられるように頑張ってください

ちなみに
自分は平均並み以下なので
デバックできていません


6.25秒以内にクリアできた方
ぜひコメントしてくださいね!!!

それ以外でも
自分のクリアタイムを
コメントしていただけると励みになります!!


今回はJavascriptを使った
簡単なゲームのご紹介でした

ぜひ遊んでみてください

それでは。





今回はベンフォードの法則と
カイ二乗検定による不正の検出方法についてです。

解説動画はこちら



ベンフォードの法則

今回はベンフォードの法則のお話です

この法則は、自然界に出てくる多くの数値の最初の桁の分布が
「一様ではなくある特定の分布になっている」という法則のことです。

電気料金の請求書、住所の番地、株価、人口の数値、死亡率、川の長さなど...
特定の範囲に限定されたものは当てはまらないですが
これに当てはまるデータが世の中にはたくさんあります。

この法則によると
大きな数値ほど最初の桁に現れる確率は小さくなり
ベンフォードの法則に従った最初の桁の理論的確率は
Pythonコードでは以下の式で表せます

math.log10((d+1)/d)



ここから先のコードをGoogle Colabで試す場合は
以下をインストールしてください
pip install japanize_matplotlib


import math
import numpy as np
import scipy.stats as stats
import japanize_matplotlib
import matplotlib.pyplot as plt

# ベンフォードの理論的な確率分布
def benford_distribution():
    return [math.log10(1 + 1 / d) for d in range(1, 10)]

benford_probs = benford_distribution()

# 棒グラフの描画
labels = [1, 2, 3, 4, 5, 6, 7, 8, 9]
bar_width = 0.35
x = np.arange(len(labels))
plt.bar(x - bar_width/2, benford_probs, width=bar_width, label='理論的確率', color='blue')
plt.xlabel('最初の桁')
plt.ylabel('確率')
plt.title('ベンフォードの法則の比較')
plt.xticks(x, labels)
plt.legend()
plt.tight_layout()
plt.show()
download


ベンフォードの法則は何に使えるのか

次のような不正の調査などで利用されているらしいです
会計データ
選挙データ
化学データ
経済データ

ここから先はどのように不正を発見するのかを
見てみましょう。



ベンフォードの法則に従うかを調査する方法

・カイ二乗検定を用いる方法

カイ二乗検定は観察されたデータと期待される
データの間の差を評価するための統計的手法のことです

データの分布が特定の比率に従っているかどうかを
検定するために使用されます

検定統計量(カイ二乗統計量)を求め、カイ二乗分布を用いて
統計量と自由度からp値(確率)を求めます

P値が一定の水準以下の場合(確率が低い)
あり得ないことが起きているとし
観察された比率が期待される比率と異なると
結論づける事ができます。


カイ二乗検定を用いたベンフォードの法則比率との比較


ここでは理論値と
観測値として、操作した比率のデータで検証していきます。

以下のコードで描画させる事ができます。
import numpy as np
import scipy.stats as stats
import japanize_matplotlib
import matplotlib.pyplot as plt

# サンプルサイズを設定
n_samples = 1000  # 適宜変更

# ベンフォードの法則に従った最初の桁の理論的確率
benford_probs = [0.301, 0.176, 0.125, 0.097, 0.079, 0.067, 0.058, 0.051, 0.046]

# 観測された出現確率を設定
observed_probs = [0.295, 0.172, 0.166, 0.092, 0.074, 0.062, 0.053, 0.046, 0.040]

# 最初の桁のラベル
labels = [1, 2, 3, 4, 5, 6, 7, 8, 9]

# 棒グラフの幅
bar_width = 0.35
x = np.arange(len(labels))

# 棒グラフの描画
plt.bar(x - bar_width/2, benford_probs, width=bar_width, label='理論的確率', color='blue')
plt.bar(x + bar_width/2, observed_probs, width=bar_width, label='観測された確率', color='orange')

# グラフの設定
plt.xlabel('最初の桁')
plt.ylabel('確率')
plt.title('ベンフォードの法則の比較')
plt.xticks(x, labels)
plt.legend()
plt.tight_layout()
plt.show()
download-1

数字の3が出てくる確率を多くしました。
これが理論上の比率と合うのかを検定します。


# 出現数を計算
observed_counts = [int(p * n_samples) for p in observed_probs]
expected_counts = [int(p * n_samples) for p in benford_probs]

print("観測数 : ",observed_counts)
print("期待値 : ",expected_counts)

# カイ二乗検定を実施
chi2_stat, p_value = stats.chisquare(f_obs=observed_counts, f_exp=expected_counts)

# 結果を表示
print("カイ二乗統計量:", chi2_stat)
print("p値:", p_value)

# 結論
alpha = 0.05  # 有意水準
if p_value < alpha:
    print("ベンフォードの法則に従わないと判断される")
else:
    print("ベンフォードの法則に従うと判断される")
観測数 :  [295, 172, 166, 92, 74, 62, 53, 46, 40]
期待値 :  [301, 176, 125, 97, 79, 67, 58, 51, 46]
カイ二乗統計量: 16.30967165997854
p値: 0.03815623468852384
ベンフォードの法則に従わないと判断される


この操作された比率の場合は
ベンフォードの法則に従わないと判断されるようです。




カイ二乗検定の仕組み



ここからはカイ二乗検定の仕組みを
説明していきます。


まずは前提として何を検定するのかを定めます。
通常は仮説を置く、ということをしています。

仮説の設定:
帰無仮説 : 観察された比率は期待される比率と等しい
対立仮説 : 観察された比率は期待される比率と異なる

理論上の比率と観測された比率が正しいのか
異なるのかを仮説に置き
一般的には帰無仮説 には差がない(比率と等しい)とします。


1.カイ二乗値を求める

仮説を置いたら
まず初めにカイ二乗値を求めていきます。

カイ二乗値は下記の式で求められます。

(観測結果 - 理論値)**2 / 理論値 の合計値

観測結果と理論値の差が大きいと
大きくなる値になります。


2.自由度を求める

次に自由度ですが、この場合は
カテゴリ数 - 1
という値になります。

この場合、1-9までの数値の比率の個数だとすると
9-1 = 8になります。


3.カイ二乗の検定統計量と比べて大きいか小さいかを見る

自由度と有意水準から
カイ二乗分布を用いて、カイ二乗検定統計量の
閾値を求めていきます。

この閾値と先ほど求めたカイ二乗値を比較します。
その結果

閾値より小さい : よく有ること
閾値より大きい : よく有ることではない

となり、閾値より大きくなった場合は
P値(確率)が低くなり

あり得ない事が起きているとして帰無仮説を棄却し
対立仮説 : 観察された比率は期待される比率と異なる
を採択することになります。


カイ二乗分布の検定統計量の表を作ってみましょう
import numpy as np
import pandas as pd
from scipy.stats import chi2

# 有意水準
alpha_levels = [0.99, 0.975, 0.95, 0.9, 0.1, 0.05, 0.025, 0.01]
# 自由度
degrees_of_freedom = range(1, 10)

# カイ二乗分布の臨界値を計算
chi_square_table = {'自由度/有意水準': degrees_of_freedom}

for alpha in alpha_levels:
    chi_square_table[f'α = {alpha}'] = [chi2.ppf(1 - alpha, df) for df in degrees_of_freedom]

# カイ二乗値のDataFrameを作成
chi_square_df = pd.DataFrame(chi_square_table)
chi_square_df
スクリーンショット 2025-01-11 17.37.20

こんな感じで、有意水準と自由度で
閾値を求められます。

紙で行う場合は、この表を見て行いますが
Pythonプログラミングでは、これを使わず
直接確率を求められます。

import numpy as np
import matplotlib.pyplot as plt
from scipy.stats import chi2

# 自由度
df = 8
# 有意水準
alpha = 0.05
# カイ二乗の臨界値を計算
critical_value = chi2.ppf(1 - alpha, df)

# 観測したカイ二乗値
observed_chi_squared = sum([(o*1000-b*1000)**2/(b*1000) for b,o in zip(benford_probs,observed_probs)])

# xの範囲を設定
x = np.linspace(0, 30, 1000)
# カイ二乗分布の確率密度関数を計算
y = chi2.pdf(x, df)

# グラフを描画
plt.figure(figsize=(10, 6))
plt.plot(x, y, label=f'Chi-squared Distribution (df={df})', color='blue')
plt.fill_between(x, y, where=(x >= critical_value), color='red', alpha=0.5, label='Rejection Region (p < 0.05)')
plt.axvline(critical_value, color='red', linestyle='--', label=f'Critical Value = {critical_value:.2f}')
plt.axvline(observed_chi_squared, color='green', linestyle='--', label=f'Observed Chi-squared = {observed_chi_squared:.2f}')
plt.title('Chi-squared Distribution with 5% Significance Level')
plt.xlabel('Chi-squared Value')
plt.ylabel('Probability Density')
plt.legend()
plt.grid()
plt.show()
download-2


この場合は求めたカイ二乗値は
閾値を超える値になるため
P値は低くなり、対立仮説の採択、となります。




実際のデータで比較するとどうなるか?

統計データを読み込みしてやってみましょう
e-statのデータを読み込みしてみます
(令和2年 都道府県・市区町村別の主な結果)

データのリンク先

こちらをColabで行う場合は
ファイル置き場においてください。


データを読み込むコードはこれです。
import pandas as pd
import numpy as np

file_path = "major_results_2020.xlsx"
  # 8行目をヘッダーとして指定(0-indexed)
df = pd.read_excel(file_path, header=8)

df.iloc[0:2, 35:]

これを集計して最初の数値のカウントと
比率を求めます。
data = df.iloc[:, 35:]

# 先頭の数字をカウントするための関数
def count_leading_digits(series):
    leading_digits = series.astype(str).str[0]
    return leading_digits.value_counts()

# 全列のデータを結合して先頭の数字をカウント
all_leading_digits = data.astype(str).stack().str[0]
digit_counts = all_leading_digits.value_counts().sort_index()
digit_counts = digit_counts[~digit_counts.index.str.contains('-')]

# 結果を表示
print(digit_counts)

# 比率を計算
total_count = digit_counts.sum()
leading_digit_ratios = digit_counts / total_count
print(leading_digit_ratios)
1    8324
2    4726
3    3465
4    2592
5    2240
6    1891
7    1594
8    1302
9    1330
Name: count, dtype: int64
1    0.303088
2    0.172080
3    0.126165
4    0.094378
5    0.081561
6    0.068854
7    0.058040
8    0.047408
9    0.048427
Name: count, dtype: float64


import numpy as np
import math
import matplotlib.pyplot as plt
from scipy.stats import norm
from collections import Counter

# ベンフォードの法則による理論的な確率分布
def benford_distribution():
    return [math.log10(1 + 1 / d) for d in range(1, 10)]

# ベンフォードの理論的分布
benford_dist = benford_distribution()

# サンプルデータの分布
sample_distribution = list(leading_digit_ratios.values)

# グラフ描画
x = range(1, 10)
fig, ax = plt.subplots(figsize=(10, 6))
ax.bar(x, sample_distribution, width=0.4, label="Sample Distribution", align="center", alpha=0.7)
ax.bar([xi + 0.4 for xi in x], benford_dist, width=0.4, label="Benford Distribution", color="orange", alpha=0.7)

# グラフの設定
ax.set_xticks(x)
ax.set_xticklabels([str(d) for d in x])
ax.set_xlabel("Leading Digit")
ax.set_ylabel("Frequency")
ax.set_title("Comparison of Sample Distribution and Benford's Law")
ax.legend()
plt.show()
download-3


カイ二乗検定を行うと
# カイ二乗検定の実施
import scipy.stats as stats

# サンプルの頻度を期待値に基づいて計算
expected_frequencies = [sum(sample_distribution) * p for p in benford_dist]
observed_frequencies = sample_distribution

# カイ二乗検定を実行
chi2_stat, p_value = stats.chisquare(observed_frequencies, f_exp=expected_frequencies)

# 結果を表示
print("カイ二乗統計量:", chi2_stat)
print("p値:", p_value)

# 有意水準を設定
alpha = 0.05
if p_value < alpha:
    print("帰無仮説を棄却します。サンプルデータはベンフォードの法則に従っていない可能性があります。")
else:
    print("帰無仮説を棄却できません。サンプルデータはベンフォードの法則に従っている可能性があります。")
カイ二乗統計量: 0.000739463432140538
p値: 0.9999999999999992
帰無仮説を棄却できません。
サンプルデータはベンフォードの法則に従っている可能性があります。


このデータの場合は
ベンフォードの法則に従っていないとは言えないようですね

正しいデータである可能性が高いです。



まとめ

ベンフォードの法則を用いると無作為に抽出したデータの場合
操作された数値を発見できる可能性はあります。

ただし、必ずしも当てはまる訳ではないので
見極めは重要です。

もしお手持ちに
会計データなどがある場合は
これで比率が正しいかを試してみると
面白いかもしれませんね。

今日はこれまでです
それでは



2つのバケツを使って水の量を測るっていう
古典的なパズル問題ありますよね

今回はそれをプログラムで解いてみました。

解説動画はこちら



二つのバケツで指定の水量を測るやつ

問題

5リットル入るバケツと3リットル入るバケツがある
この2つだけを使って4リットルの水を測りたい
どうすれば良いでしょうか

ヒント : 水は捨てたり移動したりして良い





プログラムでの解法

幅優先探索(breadth-first search, BFS)を
用いて手順の探索を行っていきます。

バケツをA, B、測りたい水量をXとし
A,Bの水量を記録していきます

途中、以下の6つの状態でうまく切り替えて手順を記録し
どこかが水量Xになったら終了になります。

1. Aを満たす
2. Bを満たす
3. Aを空にする
4. Bを空にする
5. AからBに移す
6. BからAに移す

プログラム上ではキューというデータ構造のものを使って
手順を記録していく形になります。


手順を解くコード


手順を算出するコードは次のようになります。

from collections import deque
from math import gcd

# 幅優先探索(breadth-first search, BFS)を用いて手順の探索を行う
def min_steps_to_measure(A, B, X):
    # 目標の水量がAとBの最大公約数で割り切れない場合は終了
    if X > max(A, B) or X % gcd(A, B) != 0:
        return []

    # キューとSETと状態の親を初期化
    queue, visited, parent = deque(), set(), {}

    # 初期状態をキューに追加
    queue.append((0, 0))  # (容器Aの水量, 容器Bの水量)
    parent[(0, 0)] = None  # 初期状態の親はNone

    while queue:
        a, b = queue.popleft()

        # 目標の水量に達した場合
        if a == X or b == X:
            # 結果を表示
            steps = []
            while (a, b) is not None:
                steps.append((a, b))
                # 親が存在する場合のみアンパック
                if parent[(a, b)] is not None:
                    a, b = parent[(a, b)]
                else:
                    break
            steps.reverse()
            return steps

        # 訪れた状態を記録
        if (a, b) in visited:
            continue
        visited.add((a, b))

        # 次の状態をキューに追加
        # 1. 容器Aを満たす
        if (A, b) not in visited:
            queue.append((A, b))
            parent[(A, b)] = (a, b)
        # 2. 容器Bを満たす
        if (a, B) not in visited:
            queue.append((a, B))
            parent[(a, B)] = (a, b)
        # 3. 容器Aを空にする
        if (0, b) not in visited:
            queue.append((0, b))
            parent[(0, b)] = (a, b)
        # 4. 容器Bを空にする
        if (a, 0) not in visited:
            queue.append((a, 0))
            parent[(a, 0)] = (a, b)
        # 5. 容器Aから容器Bに移す
        transfer_to_B = min(a, B - b)
        if (a - transfer_to_B, b + transfer_to_B) not in visited:
            queue.append((a - transfer_to_B, b + transfer_to_B))
            parent[(a - transfer_to_B, b + transfer_to_B)] = (a, b)
        # 6. 容器Bから容器Aに移す
        transfer_to_A = min(b, A - a)
        if (a + transfer_to_A, b - transfer_to_A) not in visited:
            queue.append((a + transfer_to_A, b - transfer_to_A))
            parent[(a + transfer_to_A, b - transfer_to_A)] = (a, b)

    return []  # 目標の水量を測ることができない場合





早速実行して見てみましょう。

バケツA,Bの値と
ターゲットの水量Xを指定すると
その手順が表示されます。

# 容器の容量と目標の水量を指定
A = 5  # 容器Aの容量
B = 3  # 容器Bの容量
X = 4  # 測りたい水量

# 最小手数を計算
steps = min_steps_to_measure(A, B, X)

# 結果を表示
if steps:
    print("途中の過程:")
    for step in steps:
        print(f"容器A: {step[0]}L, 容器B: {step[1]}L")
else:
    print("目標の水量を測ることができませんでした。")

途中の過程:
容器A: 0L, 容器B: 0L
容器A: 5L, 容器B: 0L
容器A: 2L, 容器B: 3L
容器A: 2L, 容器B: 0L
容器A: 0L, 容器B: 2L
容器A: 5L, 容器B: 2L
容器A: 4L, 容器B: 3L



このように探索の結果が表示されます。

文字だけだと分かりづらいので
描画してみることにしましょう。



手順を描画するコード

ipywidgetsでスライダーを作り
stepを可変で描画できます。

A,B,Xを色々な設定に変えて
色々遊んでみてください。

import matplotlib.pyplot as plt
import matplotlib.patches as patches
import numpy as np
import ipywidgets as widgets
from IPython.display import display
%matplotlib inline

A = 5  # 容器Aの容量
B = 3  # 容器Bの容量
X = 4  # 測りたい水量

steps = min_steps_to_measure(A, B, X)

# 描画関数
def draw_graph(step):
    a_value = steps[step][0]
    b_value = steps[step][1]
    a_frame_height = A - a_value
    b_frame_height = B - b_value

    # 描画
    fig, ax = plt.subplots(figsize=(8, 4))
    ax.bar(0, a_value, width=0.4, label='Container A Volume (L)', color='blue', linewidth=2, edgecolor='black')
    rect_a = patches.Rectangle((-0.2, a_value), 0.4, a_frame_height, linewidth=2, edgecolor='black', facecolor='none')
    ax.add_patch(rect_a)
    ax.text(0, max(A,B)+1, str(a_value), ha='center', va='bottom', fontsize=12, color='green')
    ax.bar(1, b_value, width=0.4, label='Container B Volume (L)', color='red', linewidth=2, edgecolor='black')
    rect_b = patches.Rectangle((0.8, b_value), 0.4, b_frame_height, linewidth=2, edgecolor='black', facecolor='none')
    ax.add_patch(rect_b)
    ax.text(1, max(A,B)+1, str(b_value), ha='center', va='bottom', fontsize=12, color='green')

    # タイトルとラベルの設定
    ax.set_title(f'Container A : {A} and B: {B} Target Volume : {X} with Frames step : {step}')
    ax.set_xlabel('Containers')
    ax.set_ylabel('Volume (L)')
    ax.set_xticks([0, 1])
    ax.set_xticklabels(['Container A', 'Container B'])
    ax.set_ylim(0, max(A,B)+3)
    plt.show()

# スライダーの作成
step_slider = widgets.IntSlider(value=0, min=0, max=len(steps)-1, step=1, description='Step:')
interactive_plot = widgets.interactive(draw_graph, step=step_slider)
display(interactive_plot)
step06


今回は2つのバケツを使って水量を測るやつを
プログラムにしてみました。

動画ではもう少し手順のかかるやつもやっているので
興味がある場合は見てみてください。

今回はここまでです
それでは


 

今回は日本人の年収ランキングを求めてみました

解説動画はこちら



問題

年収500万円の人は、日本の年収ランキングで
何位くらいになるでしょうか?

参考 : 日本の人口 1.245億 (2023年)





年収ランキングを求める



今回は本当に正しい値を求める事は難しいので
ザックリと近似値を求めていきます。

日本人の年収は
対数正規分布
というものに近似しています。

この分布を用いて
ざっくりと計算をしていきます。

今回は日本人の
人口が12450万人
平均年収 450万円
年収中央値を 350万円
で設定して行います。



年収分布と上位何%を計算するコード

対数正規分布と
おおよその順位を計算するコードです

Google Colabなどに貼り付けて
incomeのスライドを変えてみてください。
import numpy as np
import matplotlib.pyplot as plt
import math
import ipywidgets as widgets
from IPython.display import display
%matplotlib inline

# 対数正規分布のパラメータ設定
median_income = 350  # 中央値(単位: 万円)
mean_income = 450    # 平均値(単位: 万円)
sigma = math.sqrt(2 * math.log(mean_income / median_income)) # 標準偏差を計算
mu = np.log(median_income) # 対数正規分布の μ を計算

# x軸の範囲(年収100万円 ~ 3000万円)
x = np.linspace(100, 3000, 1000)

# パーセンタイルランクを計算する関数
def calc_percentile_rank(income):
    return 1 - 0.5 * (1 + np.sign(income - median_income) * np.sqrt(1 - np.exp(-((np.log(income) - mu) ** 2) / (2 * sigma ** 2))))

# プロットを更新する関数
def update_plot(income):
    # 確率密度関数 (PDF) を計算
    pdf = (1 / (x * sigma * np.sqrt(2 * np.pi))) * np.exp(-((np.log(x) - mu) ** 2) / (2 * sigma ** 2))
    percentile_rank = calc_percentile_rank(income) # ランク計算
    total_population = 124_500_000 # 日本の人口
    rank = (percentile_rank) * total_population
    plt.clf()
    
    # 確率密度関数をプロット
    plt.figure(figsize=(12, 4))
    plt.plot(x, pdf, label="Log-normal Distribution", color="blue")
    plt.axvline(income, color="red", linestyle="--", label=f"Income = {income}") # incomeの位置に線を引く
    plt.axvline(median_income, color='orange', linestyle='--', label=f"Median: {median_income}")
    plt.axvline(mean_income, color='green', linestyle='--', label=f"Mean: {mean_income}")
    plt.fill_between(x, pdf, where=(x >= income), color="lightblue", alpha=0.5) # income以上の範囲を薄い青で塗りつぶす
    plt.title(f"nensyu {income} : Top {percentile_rank*100:.4f}% , Rank {rank:.0f}", fontsize=14)
    plt.xlabel("Income", fontsize=12)
    plt.ylabel("Density", fontsize=12)
    plt.legend()
    plt.grid(alpha=0.5)
    plt.show()

# スライダーを作成
income_slider = widgets.IntSlider(value=500, min=100, max=3000, step=1, description='Income:')
interactive_plot = widgets.interactive(update_plot, income=income_slider)
display(interactive_plot)
download


まとめ
日本人の人口が12450万人
平均年収450万円
年収中央値350の設定では

年収500万円くらいの人の年収は
約上位32.76% : 4078万位くらいになる

年収3000万円以上の人は
0.25%くらいしかいないっぽいです
(約30万人)

自分の年収を入れて
どれくらいのランクになるか
遊んでみてくださいね

今回はこれまでです
それでは




 

このページのトップヘ