Pythonでの統計解析入門と実践的な使い方17選

import numpy as np

# サンプルデータの作成
data = np.array([1, 2, 3, 4, 5, 6, 7, 8, 9, 10])

# 平均値の計算
mean = np.mean(data)
print(f"平均値: {mean}")

# 中央値の計算
median = np.median(data)
print(f"中央値: {median}")

# 標準偏差の計算
std_dev = np.std(data)
print(f"標準偏差: {std_dev}")

# 最小値と最大値の計算
min_val = np.min(data)
max_val = np.max(data)
print(f"最小値: {min_val}")
print(f"最大値: {max_val}")

平均値: 5.5
中央値: 5.5
標準偏差: 2.8722813232690143
最小値: 1
最大値: 10

import numpy as np
import matplotlib.pyplot as plt

# より多くのサンプルデータを生成
data = np.random.normal(0, 1, 1000)  # 平均0、標準偏差1の正規分布から1000個のデータを生成

# ヒストグラムの作成
plt.hist(data, bins=30, edgecolor='black')
plt.title('データの分布')
plt.xlabel('値')
plt.ylabel('頻度')
plt.show()

import pandas as pd
import numpy as np

# サンプルデータの作成
np.random.seed(0)
dates = pd.date_range('20230101', periods=100)
df = pd.DataFrame({
    'date': dates,
    'product': np.random.choice(['A', 'B', 'C'], 100),
    'sales': np.random.randint(100, 1000, 100),
    'quantity': np.random.randint(1, 50, 100)
})

# 一部のデータを欠損値に置き換え
df.loc[10:20, 'sales'] = np.nan

print("元のデータ：")
print(df.head())
print("\n統計情報：")
print(df.describe())

# 欠損値の処理
df['sales'] = df['sales'].fillna(df['sales'].mean())

# 異常値の検出と処理（例：売上が平均+3*標準偏差を超える場合）
sales_mean = df['sales'].mean()
sales_std = df['sales'].std()
df.loc[df['sales'] > sales_mean + 3*sales_std, 'sales'] = sales_mean + 3*sales_std

# データ型の変換
df['date'] = pd.to_datetime(df['date'])
df['product'] = df['product'].astype('category')

# 新しい特徴量の作成
df['revenue'] = df['sales'] * df['quantity']

print("\nクリーニング後のデータ：")
print(df.head())
print("\n統計情報：")
print(df.describe())

元のデータ：
        date product  sales  quantity
0 2023-01-01       C    345        18
1 2023-01-02       B    773        24
2 2023-01-03       A    951         2
3 2023-01-04       B    389        16
4 2023-01-05       B    941        27

統計情報：
             sales     quantity
count   89.000000  100.000000
mean   539.134831   24.220000
std    257.246757   13.546108
min    106.000000    1.000000
25%    331.000000   13.000000
50%    543.000000   24.000000
75%    743.000000   35.250000
max    983.000000   49.000000

クリーニング後のデータ：
        date product       sales  quantity      revenue
0 2023-01-01       C  345.000000        18   6210.00000
1 2023-01-02       B  773.000000        24  18552.00000
2 2023-01-03       A  951.000000         2   1902.00000
3 2023-01-04       B  389.000000        16   6224.00000
4 2023-01-05       B  941.000000        27  25407.00000

統計情報：
             sales     quantity        revenue
count  100.000000  100.000000     100.000000
mean   539.134831   24.220000   13493.203064
std    231.507676   13.546108    9658.737542
min    106.000000    1.000000     300.000000
25%    356.250000   13.000000    6199.250000
50%    543.000000   24.000000   11899.000000
75%    720.750000   35.250000   18769.500000
max   1052.468750   49.000000   43821.000000

import pandas as pd
import matplotlib.pyplot as plt
import seaborn as sns

# 前のコードで作成したdfを使用します

# 記述統計量の算出
desc_stats = df.describe()
print("記述統計量：")
print(desc_stats)

# 製品ごとの売上の箱ひげ図
plt.figure(figsize=(10, 6))
sns.boxplot(x='product', y='sales', data=df)
plt.title('製品ごとの売上分布')
plt.show()

# 売上と数量の散布図
plt.figure(figsize=(10, 6))
sns.scatterplot(x='quantity', y='sales', hue='product', data=df)
plt.title('売上と数量の関係')
plt.show()

# 日付ごとの売上推移
df.set_index('date')['sales'].plot(figsize=(10, 6))
plt.title('日付ごとの売上推移')
plt.show()

# 相関係数の算出
corr = df[['sales', 'quantity', 'revenue']].corr()
print("\n相関係数：")
print(corr)

# 相関係数のヒートマップ
plt.figure(figsize=(8, 6))
sns.heatmap(corr, annot=True, cmap='coolwarm')
plt.title('変数間の相関係数')
plt.show()

記述統計量：
             sales     quantity        revenue
count  100.000000  100.000000     100.000000
mean   539.134831   24.220000   13493.203064
std    231.507676   13.546108    9658.737542
min    106.000000    1.000000     300.000000
25%    356.250000   13.000000    6199.250000
50%    543.000000   24.000000   11899.000000
75%    720.750000   35.250000   18769.500000
max   1052.468750   49.000000   43821.000000

相関係数：
              sales  quantity    revenue
sales     1.000000  0.045066   0.725451
quantity  0.045066  1.000000   0.701293
revenue   0.725451  0.701293   1.000000

import numpy as np
import matplotlib.pyplot as plt

# 正規分布からの乱数生成
mu, sigma = 0, 0.1  # 平均と標準偏差
normal_samples = np.random.normal(mu, sigma, 10000)

# ポアソン分布からの乱数生成
lam = 5  # 平均到着率
poisson_samples = np.random.poisson(lam, 10000)

# ヒストグラムの描画
plt.figure(figsize=(12, 5))

plt.subplot(1, 2, 1)
plt.hist(normal_samples, bins=50, density=True)
plt.title('正規分布')
plt.xlabel('値')
plt.ylabel('頻度')

plt.subplot(1, 2, 2)
plt.hist(poisson_samples, bins=range(0, 20), density=True)
plt.title('ポアソン分布')
plt.xlabel('値')
plt.ylabel('頻度')

plt.tight_layout()
plt.show()

# 統計量の計算
print("正規分布の統計量:")
print(f"平均: {np.mean(normal_samples):.4f}")
print(f"標準偏差: {np.std(normal_samples):.4f}")

print("\nポアソン分布の統計量:")
print(f"平均: {np.mean(poisson_samples):.4f}")
print(f"分散: {np.var(poisson_samples):.4f}")

正規分布の統計量:
平均: -0.0001
標準偏差: 0.1001

ポアソン分布の統計量:
平均: 4.9933
分散: 4.9859

import numpy as np
import matplotlib.pyplot as plt

# データの生成
np.random.seed(0)
X = np.random.rand(100, 1)
y = 2 + 3 * X + np.random.randn(100, 1) * 0.1

# デザイン行列の作成
X_b = np.c_[np.ones((100, 1)), X]  # バイアス項を追加

# 最小二乗法による係数の推定
theta = np.linalg.inv(X_b.T.dot(X_b)).dot(X_b.T).dot(y)

# 予測
X_new = np.array([[0], [1]])
X_new_b = np.c_[np.ones((2, 1)), X_new]
y_pred = X_new_b.dot(theta)

# 結果のプロット
plt.scatter(X, y)
plt.plot(X_new, y_pred, "r-")
plt.xlabel("X")
plt.ylabel("y")
plt.title("線形回帰モデル")
plt.show()

print("推定された係数:")
print(f"切片: {theta[0][0]:.4f}")
print(f"傾き: {theta[1][0]:.4f}")

推定された係数:
切片: 2.0024
傾き: 2.9959

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

# データの生成
np.random.seed(0)
group1 = np.random.normal(loc=50, scale=10, size=100)  # 平均50、標準偏差10の正規分布
group2 = np.random.normal(loc=55, scale=10, size=100)  # 平均55、標準偏差10の正規分布

# t検定の実行
t_statistic, p_value = stats.ttest_ind(group1, group2)

# 結果の表示
print(f"t統計量: {t_statistic:.4f}")
print(f"p値: {p_value:.4f}")

# データの可視化
plt.figure(figsize=(10, 6))
plt.hist(group1, bins=20, alpha=0.5, label='グループ1')
plt.hist(group2, bins=20, alpha=0.5, label='グループ2')
plt.xlabel('値')
plt.ylabel('頻度')
plt.title('二つのグループの分布')
plt.legend()
plt.show()

t統計量: -3.7472
p値: 0.0002

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

# データの生成
np.random.seed(0)
x = np.linspace(0, 10, 100)
y = 2 * x + 1 + np.random.normal(0, 1, 100)

# 回帰分析の実行
slope, intercept, r_value, p_value, std_err = stats.linregress(x, y)

# 結果の表示
print(f"傾き: {slope:.4f}")
print(f"切片: {intercept:.4f}")
print(f"決定係数 (R^2): {r_value**2:.4f}")
print(f"p値: {p_value:.4f}")

# データと回帰直線の可視化
plt.figure(figsize=(10, 6))
plt.scatter(x, y, label='データ')
plt.plot(x, slope * x + intercept, color='r', label='回帰直線')
plt.xlabel('x')
plt.ylabel('y')
plt.title('回帰分析')
plt.legend()
plt.show()

傾き: 2.0148
切片: 0.9668
決定係数 (R^2): 0.9349
p値: 0.0000

from sklearn import datasets
from sklearn.model_selection import train_test_split
from sklearn.preprocessing import StandardScaler
from sklearn.svm import SVC
from sklearn.metrics import classification_report, confusion_matrix
import matplotlib.pyplot as plt
import seaborn as sns

# データのロードと分割
iris = datasets.load_iris()
X = iris.data
y = iris.target
X_train, X_test, y_train, y_test = train_test_split(X, y, test_size=0.3, random_state=42)

# データの標準化
scaler = StandardScaler()
X_train = scaler.fit_transform(X_train)
X_test = scaler.transform(X_test)

# モデルの訓練
svm = SVC(kernel='rbf', random_state=42)
svm.fit(X_train, y_train)

# 予測
y_pred = svm.predict(X_test)

# 結果の評価
print(classification_report(y_test, y_pred))

# 混同行列の可視化
cm = confusion_matrix(y_test, y_pred)
plt.figure(figsize=(10,7))
sns.heatmap(cm, annot=True, fmt='d')
plt.title('Confusion Matrix')
plt.ylabel('True label')
plt.xlabel('Predicted label')
plt.show()

              precision    recall  f1-score   support

           0       1.00      1.00      1.00        16
           1       0.94      1.00      0.97        15
           2       1.00      0.93      0.97        14

    accuracy                           0.98        45
   macro avg       0.98      0.98      0.98        45
weighted avg       0.98      0.98      0.98        45

from sklearn.cluster import KMeans
from sklearn.datasets import make_blobs
import matplotlib.pyplot as plt

# データの生成
X, _ = make_blobs(n_samples=300, centers=4, cluster_std=0.60, random_state=0)

# K-meansモデルの訓練
kmeans = KMeans(n_clusters=4)
kmeans.fit(X)

# 結果の可視化
plt.figure(figsize=(10, 7))
plt.scatter(X[:, 0], X[:, 1], c=kmeans.labels_, cmap='viridis')
plt.scatter(kmeans.cluster_centers_[:, 0], kmeans.cluster_centers_[:, 1], 
            marker='x', s=200, linewidths=3, color='r')
plt.title('K-means Clustering')
plt.show()

# クラスタの中心座標
print("クラスタの中心座標:")
print(kmeans.cluster_centers_)

# 各データポイントのクラスタラベル
print("\n最初の10個のデータポイントのクラスタラベル:")
print(kmeans.labels_[:10])

クラスタの中心座標:
[[ 3.02720263  2.10955276]
 [-1.82277617 -3.04672721]
 [-3.52269546  2.65780674]
 [ 2.44429694 -2.14364446]]

最初の10個のデータポイントのクラスタラベル:
[2 2 2 2 2 2 2 2 2 2]

import pandas as pd
import numpy as np
import matplotlib.pyplot as plt
from statsmodels.tsa.arima.model import ARIMA
from statsmodels.graphics.tsaplots import plot_acf, plot_pacf

# データの生成（仮想的な月次売上データ）
np.random.seed(0)
date_rng = pd.date_range(start='2018-01-01', end='2022-12-31', freq='M')
sales = pd.Series(np.random.normal(10000, 1000, len(date_rng)) + np.arange(len(date_rng)) * 50, index=date_rng)

# トレンドと季節性を追加
sales = sales + np.sin(np.arange(len(sales)) * 2 * np.pi / 12) * 1000

# データの可視化
plt.figure(figsize=(12, 6))
plt.plot(sales)
plt.title('月次売上データ')
plt.xlabel('日付')
plt.ylabel('売上')
plt.show()

# ACFとPACFのプロット
fig, (ax1, ax2) = plt.subplots(2, 1, figsize=(12, 8))
plot_acf(sales, ax=ax1)
plot_pacf(sales, ax=ax2)
plt.show()

# ARIMAモデルの構築と予測
model = ARIMA(sales, order=(1, 1, 1))
results = model.fit()

# 将来12ヶ月の予測
forecast = results.forecast(steps=12)

# 結果の可視化
plt.figure(figsize=(12, 6))
plt.plot(sales, label='実績値')
plt.plot(forecast, label='予測値', color='red')
plt.title('売上の実績値と予測値')
plt.xlabel('日付')
plt.ylabel('売上')
plt.legend()
plt.show()

print(results.summary())

import nltk
from nltk.corpus import stopwords
from nltk.tokenize import word_tokenize
from nltk.sentiment import SentimentIntensityAnalyzer
from collections import Counter
import matplotlib.pyplot as plt

# 必要なNLTKデータのダウンロード
nltk.download('punkt')
nltk.download('stopwords')
nltk.download('vader_lexicon')

# サンプルテキスト
text = """
Python is a powerful programming language for data analysis and machine learning. 
It's easy to learn and has a large community of developers. 
Many data scientists love Python because of its simplicity and versatility.
However, some people prefer other languages like R for statistical analysis.
Overall, Python is a great choice for beginners and experts alike in the field of data science.
"""

# テキストの前処理
words = word_tokenize(text.lower())
stop_words = set(stopwords.words('english'))
filtered_words = [word for word in words if word.isalnum() and word not in stop_words]

# 頻出単語の分析
word_freq = Counter(filtered_words)
top_words = word_freq.most_common(10)

# 頻出単語のグラフ化
plt.figure(figsize=(10, 6))
plt.bar(*zip(*top_words))
plt.title('Top 10 Most Frequent Words')
plt.xlabel('Words')
plt.ylabel('Frequency')
plt.xticks(rotation=45)
plt.tight_layout()
plt.show()

# 感情分析
sia = SentimentIntensityAnalyzer()
sentiment_scores = sia.polarity_scores(text)

print("感情分析結果:")
print(f"Positive: {sentiment_scores['pos']:.2f}")
print(f"Neutral: {sentiment_scores['neu']:.2f}")
print(f"Negative: {sentiment_scores['neg']:.2f}")
print(f"Compound: {sentiment_scores['compound']:.2f}")

# 文ごとの感情分析
sentences = text.split('.')
sentence_sentiments = [(sentence.strip(), sia.polarity_scores(sentence)['compound']) for sentence in sentences if sentence.strip()]

# 文ごとの感情スコアのグラフ化
plt.figure(figsize=(12, 6))
plt.bar(range(len(sentence_sentiments)), [score for _, score in sentence_sentiments])
plt.title('Sentiment Scores for Each Sentence')
plt.xlabel('Sentence')
plt.ylabel('Compound Sentiment Score')
plt.tight_layout()
plt.show()

感情分析結果:
Positive: 0.19
Neutral: 0.72
Negative: 0.09
Compound: 0.84

import numpy as np
import matplotlib.pyplot as plt
from matplotlib.gridspec import GridSpec

# データの生成
np.random.seed(0)
x = np.linspace(0, 10, 100)
y1 = 4 + 2 * np.sin(2 * x)
y2 = 4 + 2 * np.cos(2 * x)
y3 = np.random.normal(4, 0.5, 100)

# プロットの設定
fig = plt.figure(figsize=(12, 8))
gs = GridSpec(2, 2)

# サブプロット1: 線グラフ
ax1 = fig.add_subplot(gs[0, 0])
ax1.plot(x, y1, label='sin')
ax1.plot(x, y2, label='cos')
ax1.set_title('Trigonometric Functions')
ax1.set_xlabel('x')
ax1.set_ylabel('y')
ax1.legend()

# サブプロット2: ヒストグラム
ax2 = fig.add_subplot(gs[0, 1])
ax2.hist(y3, bins=20, edgecolor='black')
ax2.set_title('Normal Distribution')
ax2.set_xlabel('Value')
ax2.set_ylabel('Frequency')

# サブプロット3: 散布図
ax3 = fig.add_subplot(gs[1, :])
scatter = ax3.scatter(y1, y2, c=x, cmap='viridis')
ax3.set_title('Scatter Plot of sin vs cos')
ax3.set_xlabel('sin(x)')
ax3.set_ylabel('cos(x)')
plt.colorbar(scatter, label='x value')

# レイアウトの調整
plt.tight_layout()
plt.show()

# 統計情報の出力
print("y1の統計情報:")
print(f"平均: {np.mean(y1):.2f}")
print(f"標準偏差: {np.std(y1):.2f}")
print(f"最小値: {np.min(y1):.2f}")
print(f"最大値: {np.max(y1):.2f}")

print("\ny2の統計情報:")
print(f"平均: {np.mean(y2):.2f}")
print(f"標準偏差: {np.std(y2):.2f}")
print(f"最小値: {np.min(y2):.2f}")
print(f"最大値: {np.max(y2):.2f}")

print("\ny3の統計情報:")
print(f"平均: {np.mean(y3):.2f}")
print(f"標準偏差: {np.std(y3):.2f}")
print(f"最小値: {np.min(y3):.2f}")
print(f"最大値: {np.max(y3):.2f}")

y1の統計情報:
平均: 4.00
標準偏差: 1.41
最小値: 2.00
最大値: 6.00

y2の統計情報:
平均: 4.00
標準偏差: 1.41
最小値: 2.00
最大値: 6.00

y3の統計情報:
平均: 3.98
標準偏差: 0.48
最小値: 2.73
最大値: 5.34

import numpy as np
import pandas as pd
import seaborn as sns
import matplotlib.pyplot as plt

# データの生成
np.random.seed(0)
n_samples = 1000
df = pd.DataFrame({
    'x': np.random.normal(0, 1, n_samples),
    'y': np.random.normal(0, 1, n_samples),
    'category': np.random.choice(['A', 'B', 'C'], n_samples),
    'value': np.random.uniform(0, 10, n_samples)
})

# プロットの設定
sns.set(style="whitegrid")
fig, axes = plt.subplots(2, 2, figsize=(15, 15))

# サブプロット1: 散布図とカーネル密度推定
sns.scatterplot(data=df, x='x', y='y', hue='category', ax=axes[0, 0])
sns.kdeplot(data=df, x='x', y='y', levels=5, color="r", linewidths=1, ax=axes[0, 0])
axes[0, 0].set_title('Scatter Plot with KDE')

# サブプロット2: ボックスプロット
sns.boxplot(data=df, x='category', y='value', ax=axes[0, 1])
axes[0, 1].set_title('Box Plot')

# サブプロット3: バイオリンプロット
sns.violinplot(data=df, x='category', y='value', ax=axes[1, 0])
axes[1, 0].set_title('Violin Plot')

# サブプロット4: ペアプロット
sns.pairplot(df, vars=['x', 'y', 'value'], hue='category', height=3)
plt.tight_layout()
plt.show()

# 統計情報の出力
print(df.groupby('category').describe())

import pandas as pd
import numpy as np

# サンプルデータの作成
df = pd.DataFrame({
    'A': [1, 2, '3', 4, 5],
    'B': [10, 20, 30, '40', 50]
})

print("元のデータフレーム:")
print(df)

# エラーが発生するケース
try:
    result = df['A'].mean()
    print("列Aの平均:", result)
except TypeError as e:
    print("エラー発生:", str(e))

# 対処法1: データ型の変換
df['A'] = pd.to_numeric(df['A'], errors='coerce')
df['B'] = pd.to_numeric(df['B'], errors='coerce')

print("\nデータ型変換後のデータフレーム:")
print(df)

# 対処法2: 数値のみを抽出して計算
numeric_A = df['A'].dropna()
numeric_B = df['B'].dropna()

print("\n列Aの平均:", numeric_A.mean())
print("列Bの平均:", numeric_B.mean())

# データ型の確認
print("\nデータ型:")
print(df.dtypes)

元のデータフレーム:
   A   B
0  1  10
1  2  20
2  3  30
3  4  40
4  5  50

エラー発生: could not convert string to float: '3'

データ型変換後のデータフレーム:
     A     B
0  1.0  10.0
1  2.0  20.0
2  3.0  30.0
3  4.0  40.0
4  5.0  50.0

列Aの平均: 3.0
列Bの平均: 30.0

データ型:
A    float64
B    float64
dtype: object

import pandas as pd
import numpy as np
import matplotlib.pyplot as plt

# サンプルデータの作成
np.random.seed(0)
df = pd.DataFrame({
    'A': np.random.rand(1000),
    'B': np.random.rand(1000)
})

# 欠損値の追加
df.loc[np.random.choice(df.index, 100, replace=False), 'A'] = np.nan

print("欠損値の数:")
print(df.isnull().sum())

# 方法1: 欠損値の削除
df_dropna = df.dropna()

# 方法2: 平均値での補完
df_mean = df.fillna(df.mean())

# 方法3: 線形補間
df_interpolate = df.interpolate()

# 結果の可視化
fig, axes = plt.subplots(2, 2, figsize=(15, 10))

axes[0, 0].scatter(df['A'], df['B'])
axes[0, 0].set_title('元のデータ')

axes[0, 1].scatter(df_dropna['A'], df_dropna['B'])
axes[0, 1].set_title('欠損値削除後')

axes[1, 0].scatter(df_mean['A'], df_mean['B'])
axes[1, 0].set_title('平均値補完後')

axes[1, 1].scatter(df_interpolate['A'], df_interpolate['B'])
axes[1, 1].set_title('線形補間後')

plt.tight_layout()
plt.show()

# 各方法の統計量比較
methods = ['元データ', '欠損値削除', '平均値補完', '線形補間']
for method, data in zip(methods, [df, df_dropna, df_mean, df_interpolate]):
    print(f"\n{method}の統計量:")
    print(data['A'].describe())

import pandas as pd
import numpy as np
import os

# 大規模データセットの生成（注意：実行にはかなりの時間とディスク容量が必要です）
def create_large_csv(filename, size_mb):
    chunk_size = 100000
    chunks = []
    total_rows = size_mb * 1024 * 1024 // 100  # 1行あたり約100バイトと仮定

    for i in range(0, total_rows, chunk_size):
        chunk = pd.DataFrame({
            'id': range(i, min(i+chunk_size, total_rows)),
            'value': np.random.rand(min(chunk_size, total_rows-i))
        })
        chunks.append(chunk)
        if len(chunks) == 10:
            pd.concat(chunks).to_csv(filename, mode='a', index=False, header=False)
            chunks = []

    if chunks:
        pd.concat(chunks).to_csv(filename, mode='a', index=False, header=False)

# 1GB程度のCSVファイルを作成
filename = 'large_dataset.csv'
if not os.path.exists(filename):
    create_large_csv(filename, 1000)  # 1000MB = 1GB

# メモリエラーが発生する可能性がある処理
try:
    df = pd.read_csv(filename)
    print("全データの平均:", df['value'].mean())
except MemoryError:
    print("メモリエラーが発生しました。チャンク処理を使用します。")

# チャンク処理による対処法
chunk_size = 1000000  # 100万行ずつ処理
total_sum = 0
total_count = 0

for chunk in pd.read_csv(filename, chunksize=chunk_size):
    total_sum += chunk['value'].sum()
    total_count += len(chunk)

average = total_sum / total_count
print("チャンク処理による全データの平均:", average)

# データタイプの最適化
optimized_dtypes = {'id': 'int32', 'value': 'float32'}
total_sum = 0
total_count = 0

for chunk in pd.read_csv(filename, chunksize=chunk_size, dtype=optimized_dtypes):
    total_sum += chunk['value'].sum()
    total_count += len(chunk)

optimized_average = total_sum / total_count
print("データタイプ最適化後の全データの平均:", optimized_average)

# メモリ使用量の比較
import psutil

def memory_usage():
    process = psutil.Process(os.getpid())
    return process.memory_info().rss / 1024 ** 2  # MBで返す

print(f"プログラム終了時のメモリ使用量: {memory_usage():.2f} MB")

import pandas as pd
import numpy as np
from sklearn.preprocessing import StandardScaler
from sklearn.cluster import KMeans
import matplotlib.pyplot as plt
import seaborn as sns

# サンプルデータの生成
np.random.seed(0)
n_customers = 1000
last_purchase = pd.to_datetime('2023-06-01') - pd.to_timedelta(np.random.randint(1, 365, n_customers), unit='D')
frequency = np.random.randint(1, 20, n_customers)
monetary = np.random.randint(10, 1000, n_customers)

df = pd.DataFrame({
    'customer_id': range(n_customers),
    'last_purchase': last_purchase,
    'frequency': frequency,
    'monetary': monetary
})

# RFM指標の計算
df['recency'] = (pd.to_datetime('2023-06-01') - df['last_purchase']).dt.days
df['frequency'] = df['frequency']
df['monetary'] = df['monetary']

# データの標準化
scaler = StandardScaler()
df_scaled = scaler.fit_transform(df[['recency', 'frequency', 'monetary']])

# K-meansクラスタリング
kmeans = KMeans(n_clusters=4, random_state=42)
df['cluster'] = kmeans.fit_predict(df_scaled)

# 結果の可視化
fig, axes = plt.subplots(2, 2, figsize=(15, 15))

sns.scatterplot(data=df, x='recency', y='frequency', hue='cluster', ax=axes[0, 0])
axes[0, 0].set_title('Recency vs Frequency')

sns.scatterplot(data=df, x='recency', y='monetary', hue='cluster', ax=axes[0, 1])
axes[0, 1].set_title('Recency vs Monetary')

sns.scatterplot(data=df, x='frequency', y='monetary', hue='cluster', ax=axes[1, 0])
axes[1, 0].set_title('Frequency vs Monetary')

sns.boxplot(data=df, x='cluster', y='monetary', ax=axes[1, 1])
axes[1, 1].set_title('Monetary by Cluster')

plt.tight_layout()
plt.show()

# クラスタごとの統計情報
print(df.groupby('cluster')[['recency', 'frequency', 'monetary']].mean())

        recency  frequency   monetary
cluster                              
0       183.682     10.122  508.39024
1        88.809     10.191  504.31915
2       270.868     10.075  503.71698
3        10.325      9.992  502.25397

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

# データの生成
np.random.seed(42)
n_A = 1000  # グループAのサンプルサイズ
n_B = 1000  # グループBのサンプルサイズ
conv_rate_A = 0.10  # グループAの基準変換率
conv_rate_B = 0.12  # グループBの変換率（効果があると仮定）

# シミュレーションデータの生成
conversions_A = np.random.binomial(1, conv_rate_A, n_A)
conversions_B = np.random.binomial(1, conv_rate_B, n_B)

# 変換率の計算
obs_conv_rate_A = np.mean(conversions_A)
obs_conv_rate_B = np.mean(conversions_B)

print(f"観測された変換率 A: {obs_conv_rate_A:.4f}")
print(f"観測された変換率 B: {obs_conv_rate_B:.4f}")

# 統計的検定（カイ二乗検定）
contingency_table = np.array([[sum(conversions_A), n_A - sum(conversions_A)],
                              [sum(conversions_B), n_B - sum(conversions_B)]])
chi2, p_value = stats.chi2_contingency(contingency_table)[:2]

print(f"カイ二乗統計量: {chi2:.4f}")
print(f"p値: {p_value:.4f}")

# 効果量（相対的な変化）の計算
relative_change = (obs_conv_rate_B - obs_conv_rate_A) / obs_conv_rate_A
print(f"相対的な変化: {relative_change:.2%}")

# 信頼区間の計算
z = stats.norm.ppf(0.975)  # 95%信頼区間のz値
se = np.sqrt(obs_conv_rate_A*(1-obs_conv_rate_A)/n_A + obs_conv_rate_B*(1-obs_conv_rate_B)/n_B)
margin_of_error = z * se
ci_lower = relative_change - margin_of_error
ci_upper = relative_change + margin_of_error

print(f"95%信頼区間: [{ci_lower:.2%}, {ci_upper:.2%}]")

# 結果の可視化
plt.figure(figsize=(10, 6))
plt.bar(['A', 'B'], [obs_conv_rate_A, obs_conv_rate_B], yerr=[se, se], capsize=10)
plt.title('A/Bテスト結果')
plt.ylabel('変換率')
plt.ylim(0, max(obs_conv_rate_A, obs_conv_rate_B) * 1.2)
plt.text('A', obs_conv_rate_A, f'{obs_conv_rate_A:.2%}', ha='center', va='bottom')
plt.text('B', obs_conv_rate_B, f'{obs_conv_rate_B:.2%}', ha='center', va='bottom')
plt.show()

観測された変換率 A: 0.0970
観測された変換率 B: 0.1150
カイ二乗統計量: 1.8121
p値: 0.1782
相対的な変化: 18.56%
95%信頼区間: [-7.74%, 44.85%]

import pandas as pd
import numpy as np
from statsmodels.tsa.statespace.sarimax import SARIMAX
import matplotlib.pyplot as plt

# データの生成
np.random.seed(0)
date_rng = pd.date_range(start='2018-01-01', end='2023-05-31', freq='M')
n = len(date_rng)
trend = np.linspace(1000, 1500, n)
seasonal = 200 * np.sin(np.arange(n) * 2 * np.pi / 12)
noise = np.random.normal(0, 50, n)
sales = trend + seasonal + noise

df = pd.DataFrame(data={'date': date_rng, 'sales': sales})
df.set_index('date', inplace=True)

# モデルの構築
model = SARIMAX(df['sales'], order=(1, 1, 1), seasonal_order=(1, 1, 1, 12))
results = model.fit()

# 将来の予測
forecast = results.get_forecast(steps=12)
forecast_ci = forecast.conf_int()

# 結果の可視化
plt.figure(figsize=(12, 6))
plt.plot(df.index, df['sales'], label='実績値')
plt.plot(forecast.predicted_mean.index, forecast.predicted_mean, color='r', label='予測値')
plt.fill_between(forecast_ci.index, forecast_ci.iloc[:, 0], forecast_ci.iloc[:, 1], color='r', alpha=0.1)
plt.title('月次販売データの需要予測')
plt.xlabel('日付')
plt.ylabel('販売量')
plt.legend()
plt.show()

# モデルの性能評価
print(results.summary())

# 予測値
print("\n今後12ヶ月の予測値:")
print(forecast.predicted_mean)

2023-06-30    1676.849741
2023-07-31    1691.925246
2023-08-31    1730.930659
...

import statsmodels.api as sm

# ARIMAモデルの例
data = sm.datasets.sunspots.load_pandas().data
y = data['SUNACTIVITY']
model = sm.tsa.ARIMA(y, order=(2,1,0))
results = model.fit()
print(results.summary())

import pandas as pd
from sklearn.model_selection import train_test_split
from sklearn.ensemble import RandomForestClassifier
from sklearn.metrics import accuracy_score

# データの読み込み
train_data = pd.read_csv('train.csv')
test_data = pd.read_csv('test.csv')

# 特徴量とターゲットの分離
features = ['Pclass', 'Sex', 'Age', 'SibSp', 'Parch', 'Fare']
X = train_data[features]
y = train_data['Survived']

# データの前処理
X['Sex'] = X['Sex'].map({'male': 0, 'female': 1})
X['Age'].fillna(X['Age'].median(), inplace=True)
X['Fare'].fillna(X['Fare'].median(), inplace=True)

# モデルのトレーニングと評価
X_train, X_val, y_train, y_val = train_test_split(X, y, test_size=0.2, random_state=42)
model = RandomForestClassifier(n_estimators=100, random_state=42)
model.fit(X_train, y_train)
predictions = model.predict(X_val)
print(f"Validation Accuracy: {accuracy_score(y_val, predictions)}")

# テストデータの予測
test_data['Sex'] = test_data['Sex'].map({'male': 0, 'female': 1})
test_data['Age'].fillna(test_data['Age'].median(), inplace=True)
test_data['Fare'].fillna(test_data['Fare'].median(), inplace=True)
test_predictions = model.predict(test_data[features])