サムネがコーヒーの記事は書きかけです。

【研究用スクリプト】複数実験系の頻度分布のカーネル密度推定を3次元空間にプロットするスクリプト

実験系が複数あって、それぞれの実験系でのデータが離散的、かつ一次元の場合は以下の方法を使用してデータを一つの画像にまとめることができます。

今回は、化学物質への曝露時間に対応する細胞内の異常構造形成数を定量化する想定でグラフを作成します。

import numpy as np
import pandas as pd
from matplotlib import pyplot as plt
import seaborn as sns
from scipy.stats import gaussian_kde
# カイ二乗分布の自由度
degrees_of_freedom = 1
# ランダムデータ生成
random_data = np.random.chisquare(degrees_of_freedom, 1000)
# 0から100の範囲にスケーリング
scaled_data = random_data * 100 / np.max(random_data)
# リストに変換
data_list_Chi = scaled_data.tolist()
# ベータ分布のパラメータ
alpha = 2  # 形状パラメータ1
beta = 5   # 形状パラメータ2
# ランダムデータ生成
random_data = np.random.beta(alpha, beta, 1000)
# 0から100の範囲にスケーリング
scaled_data = random_data * 100
# リストに変換
data_list_beta = scaled_data.tolist()
# ポアソン分布の平均
lambda_param = 5
# ランダムデータ生成
random_data = np.random.poisson(lambda_param, 1000)
# 0から100の範囲にクリップ
clipped_data = np.clip(random_data, 0, 100)
# リストに変換
data_list_poisson = clipped_data.tolist()
# ガンマ分布のパラメータ
shape_param = 2.0
scale_param = 10.0
# ランダムデータ生成
random_data = np.random.gamma(shape_param, scale=scale_param, size=1000)
# 0から100の範囲にクリップ
clipped_data = np.clip(random_data, 0, 100)
# リストに変換
data_list_gamma = clipped_data.tolist()
df = pd.DataFrame({
    '2.0': [np.random.randint(0,100) for i in range(1000)],
    '1.5': data_list_Chi,
    '1.0': data_list_beta,
    '0.5': [np.random.random()*100 for i in range(1000)],
    '0.0': data_list_gamma
})
# カーネル密度推定の範囲を設定 (0以上)
density_range = np.linspace(0, 150, 500)  # 適切な範囲を指定
# 3Dプロットの設定
fig = plt.figure()
ax = fig.add_subplot(111, projection='3d')
# 各カラムのヒストグラムとカーネル密度推定を描画
for j, column in enumerate(df.columns):
    data = df[column]
    data = data[data >= 0]  # 0以上のデータに絞る
    density = gaussian_kde(data)
    
    # 密度関数の計算
    density_values = density(density_range)
    
    # 端を繋げるための補完
    density_values[:10] = density_values[10]  # 先頭の値をコピー
    density_values[-10:] = density_values[-11]  # 末尾の値をコピー
    
    zs = np.full_like(density_range, [120,90,60,30,0][j])
    cmap = plt.get_cmap('viridis')
    colors = [cmap(i / 4) for i in range(5)]
    ax.plot(density_range, zs, density_values, color = colors[j])
# 軸ラベルの設定
ax.set_ylabel('Exposure time (min)')
ax.set_zlabel('Density')
ax.set_xlabel('Num aggragation')
fig.savefig("result.png",dpi = 500)
>>>
None

実行結果

コメントを残す

メールアドレスが公開されることはありません。 が付いている欄は必須項目です