Copying "Introduction to Multivariate Analysis" in Python Vol. 4 - Chapter 2 "Introduction to Multiple Regression Analysis" (3) Statistical Testing in Multiple Regression Analysis

Chapter 2 "Introduction to Multiple Regression Analysis"

Authors of the book: Dr. Sadao Ishimura, Dr. Koshiro Ishimura

This is a record of Python coding for Chapter 2 "Introduction to Multiple Regression Analysis" of the book "Introduction to Multivariate Analysis".

This is a series where we learn the basics of multivariate analysis together with Python.

This article covers statistical testing in multiple regression analysis. Specifically, we will work on
the ANOVA table for multiple regression (testing the significance of regression) and testing the significance of partial regression coefficients.

We will proceed with ChatGPT-assisted learning!
Now, let's open the book and embark on a journey of multivariate analysis 🚀

Introduction

This blog series introduces the "joy of multivariate analysis" gained through Python coding of the book "Introduction to Multivariate Analysis" (Tokyo Tosho, referred to as the "text").

Information on the book and citation notation are posted in the linked article.

Chapter 2: Introduction to Multiple Regression Analysis

This article covers the following section of Chapter 2.

The data used in this article is cited directly from the data published in the text.
For data with a small number of entries, we register the data in the code, and for data with a large number of entries, we create a CSV file and load it.

We will import the libraries used in this article.

### インポート

# 数値計算
import numpy as np
import pandas as pd

# 統計
import scipy.stats as stats
import statsmodels.api as sm
import statsmodels.formula.api as smf

# 描画
import matplotlib.pyplot as plt
plt.rcParams['font.family'] = 'Meiryo'  # または import japanize_matplotlib

Preparation for Analysis

We will cite Table 2.1.1 "Data on Temperature, Pressure, and Orientation" from page 35 of the text.

■ Loading Data
Set the data in a pandas DataFrame.

### 温度・圧力・配向度のデータ p.35 表2.1.1

# データの登録
data1 = pd.DataFrame(
    {'配向度': [45, 38, 41, 34, 59, 47, 35, 43, 54, 52],
    '温度': [17.5, 17.0, 18.5, 16.0, 19.0, 19.5, 16.0, 18.0, 19.0, 19.5],
    '圧力': [30, 25, 20, 30, 45, 35, 25, 35, 35, 40],
    }, index=range(1, 11))
data1.index.name = 'サンプルNo.'

# 結果の表示
data1

【Execution Results】
The number of data points (sample size) is 10.
The objective variable is "Orientation," and the candidate explanatory variables are "Temperature" and "Pressure."

■ Common Settings
Set the explanatory variable names and the objective variable name.
This is useful when identifying columns in the DataFrame.

## 共通設定

# 説明変数 
VARS = ['温度', '圧力']
# 目的変数
TARGET = '配向度'

[Execution Result] None

■ Obtaining Multiple Regression Analysis Results with Python Libraries
We will perform a multiple regression analysis using the statsmodels library, using "Temperature" and "Pressure" as explanatory variables.

🖲️statsmodels
Set the relationship between the objective variable and explanatory variables using the formula argument.
The basic syntax is "objective variable ~ explanatory variable 1 + explanatory variable 2 + ...".

## 重回帰分析の実行 statsmodels利用

result1_sm = smf.ols(formula='配向度 ~ 温度 + 圧力', data=data1).fit()
result1_sm.summary()

[Execution Result]
"result1_sm", which stores the results of the multiple regression analysis, will appear frequently in this article!
The summary of the multiple regression analysis is as follows.

Statistical Testing p.64

The text introduces the following two types of statistical tests for multiple regression analysis:
・Test for significance of regression (ANOVA table for multiple regression)
・Test for significance of partial regression coefficients (test for explanatory variables)

Statistical testing is "using sample data to determine whether a hypothesis about a population is statistically valid."
Generally, it is executed in the following steps (partially modified from the text's steps).

1. Set the null hypothesis $${H_0}$$ and alternative hypothesis $${H_1}$$ regarding the population characteristic values

2. Calculate the test statistic from the sample data

3. Reject the null hypothesis $${H_0}$$ when the test statistic is included in the rejection region

Test for significance of regression (ANOVA table for multiple regression) p.65~

This is a statistical test regarding whether the "entire multiple regression model" has an effect on the objective variable.

① Hypothesis
■ Text version:
・Null hypothesis $${H_0}$$: The multiple regression equation is not useful for prediction

■ Common version (regarding $${p}$$ population partial regression coefficients):
・Null hypothesis $${H_0}$$: $${\beta_1=\beta_2 = \cdots = \beta_p=0}$$
　(All population partial regression coefficients are 0)
・Alternative hypothesis $${H_1}$$: $${\beta_1 \neq 0\ \text{or}\ \beta_2 \neq 0\ \text{or} \cdots \text{or}\ \beta_p \neq 0}$$
　(At least one population partial regression coefficient is not 0)

② Test statistic
The test statistic used in the regression significance test is the $${F}$$ value.
The text calculates the $${F}$$ value using an "ANOVA table" from the statistical analysis software SPSS.

Let's check how to calculate the $${F}$$ value using an ANOVA table with calculation steps.

$$
\begin{array}{|c|c:c:c|c|}
\hline
\\
Source of Variation & Sum of Squares & Degrees of Freedom & Mean Square (Variance) & F-value \\
\\
\hline
\\
Regression & S_R & p & V_R = S_R / p & F = V_R / V_e \\
\\
\hdashline
\\
Residual & S_E & N-p-1 & V_e = S_E / (N-p-1) & - \\
\\
\hline
\\
Total & S_T & N-1 & - & - \\
\\
\hline
\end{array}
$$

It is calculated using the regression sum of squares $${S_R}$$, the residual sum of squares $${S_E}$$, the number of data points (sample size) $${N}$$, and the number of explanatory variables $${p}$$.
The sum of squares was introduced in the article before last.

And the $${F}$$ value follows an "$${F}$$ distribution with degrees of freedom $${(p,\ N-p-1)}$$".

🔢 Creating and calculating the ANOVA table
Let's calculate the ANOVA table using the regression analysis results from statsmodels.
Pass the regression analysis result result1_sm to the anova_lm function.

### 重回帰の分散分析表 F検定 p.65 表2.6.2 ※statsmodels利用
sm.stats.anova_lm(result1_sm, typ=1)

[Execution Result]

The regression rows are separated by explanatory variable, so I cannot check the F-value for the entire regression...

So, let's create a function to generate an ANOVA table using the regression analysis results from statsmodels! (DIY)

### 重回帰の分散分析表作成関数 for statsmodels p.65

def anova_table(result):
    ## 設定と準備
    # 標本サイズn、説明変数の数p
    n, p = result.nobs, result.df_model
    # 欠損値
    NaN = np.nan
    
    ## データフレームの作成
    # statsmodelsのOLSResultから平均和、自由度、平均平方、F値、P値を設定
    df = pd.DataFrame({
        '平方和': [result.ess, result.ssr, result.ess + result.ssr],
        '自由度': [int(p), int(n - p - 1), int(n - 1)],
        '平均平方': [result.mse_model, result.mse_resid, NaN],
        'F値': [result.fvalue, NaN, NaN],
        'p値': [result.f_pvalue, NaN, NaN]},
        index=['回帰', '残差', '全体'])
    
    ## 戻り値: 分散分析表データフレーム
    return df

Execute the function by passing the statsmodels regression analysis result result1_sm as an argument.

# 分散分析表の作成
anova_table(result1_sm).round(3)

[Execution Result]
It now looks the same as Table 2.6.2 "ANOVA table by SPSS" on page 65 of the textbook.
The F-value is 21.176.
The p-value is 0.001 <= 0.05.

Actually, you can check the F-value and its p-value by displaying the "summary" of the statsmodels regression analysis results without creating an ANOVA table.

# statsmodels のF値
result1_sm.summary().tables[0]

[Execution Result]
The F-value is F-statistic = 21.18, and the p-value is Prob (F-statistic) = 0.00107.

By the way, to extract only the F-value:

# statsmodels のF値
result1_sm.fvalue

[Execution Result]

Now, let's move on to judging the hypothesis.

(3) Judging whether the null hypothesis can be rejected
To state the conclusion first...
The p-value <= significance level of 5%.
At a 5% significance level, the null hypothesis is rejected, and we can say that the multiple regression model is significant.

A multiple regression model using "temperature" and "pressure" as explanatory variables is useful for predicting the degree of orientation.
It was worth creating the multiple regression model!!!

Let's perform a visualization corresponding to Figure 2.6.1 "Significance Probability and Significance Level" on page 65 of the textbook.
I will create a drawing function.

### F分布の可視化 for statsmodels

def plot_f_dist(result, alpha=0.05):
    
    ## 設定と準備
    # 基本色の設定
    color = 'tab:blue'
    # 標本サイズn、説明変数の数p
    n, p = result.nobs, result.df_model
    # 重回帰モデルのF値とP値の取得
    f_val, p_val = result.fvalue, result.f_pvalue
    # F分布の自由度の算出
    dfn, dfd = int(p), int(n-p-1)
    # F分布の設定
    f_dist = stats.f(dfn=dfn, dfd=dfd)
    # 棄却限界値の算出
    c_val = f_dist.isf(q=alpha)
    # 描画のためのx軸の値
    x_val1 = np.linspace(0, f_val * 1.1)
    x_val2 = np.linspace(c_val, f_val * 1.1)
    
    ## 描画
    # F分布の確率密度関数の描画
    plt.plot(x_val1, f_dist.pdf(x_val1), color=color, label='$F$分布')
    # 有意水準alphaの棄却域の塗りつぶし描画
    plt.fill_between(x_val2, 0, f_dist.pdf(x_val2), color=color, alpha=0.2,
                     label='棄却域')
    # 棄却限界値の垂直点線（黒）の描画
    plt.axvline(c_val, color='black', ls='--', label='棄却限界値')
    # F値の垂直点線の描画
    plt.axvline(f_val, color='tab:red', ls='--', label='F値')
    # 修飾
    plt.title(f'自由度{dfn, dfd}の$F$分布\n'
              f'棄却限界値:{c_val:.4f}, F値:{f_val:.4f}'
              f'\n有意水準{alpha:.1%}, $p$値={p_val:.4f}')
    plt.xticks([c_val, f_val])
    plt.legend(loc='upper right')
    plt.show()

Now, let's draw it.

### 有意水準の可視化 p.65 図2.6.1
plot_f_dist(result1_sm)

[Execution Result]
If the F-value is greater than or equal to the rejection limit, the null hypothesis is rejected.
In this chart, the F-value of 21.1756 exceeds the rejection limit of 4.7374, so we can see that the null hypothesis can be rejected.

🛸 A little detour: ChatGPT recommended chart 🛸

I wanted to know more examples of multiple regression model visualization, so I asked ChatGPT.

The answer is... 🥁

Simulation distribution of "explanatory dominance (R²)" under the null hypothesis 💹

It is a chart that allows you to grasp the effectiveness of a multiple regression model using the coefficient of determination R^2.
Let's draw it for now!

### 決定係数 R² の帰無分布をシミュレーション

# R² の帰無分布は、帰無仮説「説明変数と目的変数に何の関係もない」のもとで、
# 偶然でどのくらいの R² が出るかを集めた分布です。
# 観測データで得られた実際の R² が、帰無分布でほとんど現れないほど大きければ、
# 「偶然ではこんなに高い説明力はほぼ起こらない」
# ⇒「モデルは意味がある（有意）」と判断します。

## 設定と準備
# シミュレーション回数
n_sim = 1000 
# 重回帰分析の決定係数
observed_r2 = result1_sm.rsquared
# 乱数生成器
rng = np.random.default_rng(seed=123)

## 帰無分析のシミュレーション
# 帰無分布の決定係数を格納するリストの初期化
r2_null = []
# シミュレーションの実行
for _ in range(n_sim):
    # 目的変数をランダムに並び替える
    y_perm = rng.permutation(data1[TARGET].values)
    # 並び替えた目的変数と並び替えていない説明変数で重回帰分析を実行
    perm_model = sm.OLS(y_perm, sm.add_constant(data1[VARS])).fit()
    # 決定係数を格納
    r2_null.append(perm_model.rsquared)

## 描画
# 描画領域の設定
plt.figure(figsize=(6, 4))
# 帰無分布の決定係数のヒストグラムの描画
plt.hist(r2_null, bins=30, density=True, edgecolor='white', alpha=0.7)
plt.axvline(observed_r2, ls='--', lw=2, color='tab:red', 
            label=f'観測 $R^2$ = {observed_r2:.2f}')
plt.xlabel('帰無分布のもとの $R^2$：帰無仮説「温度・圧力と配向度は無関係」',
           fontsize=12)
plt.ylabel('密度', fontsize=12)
plt.title('$R^2$ の帰無分布シミュレーションと観測値')
plt.legend(loc='upper right')
plt.tight_layout()
plt.show()

[Execution Result]

This histogram shows the distribution for the case where the "null hypothesis: there is no relationship between the explanatory variables and the objective variable" is true.
Specifically, it repeats the trial of performing multiple regression analysis and calculating the coefficient of determination 1000 times using "randomly shuffled objective variables" and "original, unshuffled explanatory variables".
It seems unlikely that the shuffled objective variables and the original explanatory variables have any relationship, right?
Reading the histogram, we can see that the coefficient of determination for most of the data is less than $${0.8}$$.

And the coefficient of determination (observed $${R^2}$$) obtained from the actual data is $${0.86}$$, indicated by the red vertical dotted line.

The distribution for the case where there is no relationship between the explanatory variables and the objective variable shows that a coefficient of determination of $${0.86}$$ rarely occurs.
"Rarely occurs ⇒ There is a relationship between the explanatory variables and the objective variable" can be considered.

And then, ChatGPT also...
taught me the relationship between the $${F}$$ value and $${R^2}$$.

From the relationship of Regression Sum of Squares $${S_R\ +}$$ Residual Sum of Squares $${S_E\ =\ }$$ Total Sum of Squares $${S_T}$$, we get

$$
\begin{align*}
R^2 = \cfrac{S_R}{S_T} &\Longrightarrow S_R = R^2 \times S_T \\
R^2 = 1-\cfrac{S_E}{S_T} &\Longrightarrow \ S_E = (1-R^2) \times S_T
\end{align*}
$$

and by substituting $${S_R, S_T}$$ into the formula for the $${F}$$ value, we get

$$
F=\cfrac{\cfrac{S_R}{p}}{\cfrac{S_E}{N-p-1}} = \cfrac{\cfrac{R^2 S_T}{p}}{\cfrac{(1-R^2)S_T}{N-p-1}} = \cfrac{\cfrac{R^2}{p}}{\cfrac{1-R^2}{N-p-1}}
$$

which results in this.

I understood the relationship that a larger coefficient of determination $${R^2}$$ leads to a larger $${F}$$ value.

Testing the significance of partial regression coefficients (testing explanatory variables) p.66~

This is a statistical test regarding whether "individual explanatory variables" have an effect on the objective variable.

(1) Hypothesis
This is a hypothesis regarding the population partial regression coefficient $${\beta_i}$$ of a specific explanatory variable.
■ Text version:
・Null hypothesis $${H_0}$$: Explanatory variables have no effect on the objective variable

■ Common version
・Null hypothesis $${H_0}$$: $${\beta_i=0}$$
　(The population partial regression coefficient $${\beta_i}$$ is 0)
・Alternative hypothesis $${H_1}$$: $${\beta_i \neq 0}$$
　(The population partial regression coefficient $${\beta_i}$$ is not 0)

(2) Test statistic
The test statistic used in the significance test of partial regression coefficients is the $${t}$$ value.
The text calculates the $${t}$$ value using the statistical analysis software SPSS.

🔢 Calculation of test statistic $${t}$$ value
Let's check the $${t}$$ value using the regression analysis results from statsmodels.

### 検定統計量の算出 p.67 表2.6.3 ※statsmodels利用
result1_sm.summary().tables[1]

[Execution Result]
The $${t}$$ value for each partial regression coefficient is displayed in the t column.

By the way, the $${t}$$ value

can be calculated with this.

🛸 A little detour: Formula for the $${\boldsymbol{t}}$$ value 🛸

The formula for the $$${t}$$$ value was not included in the text, so I was curious...
I asked ChatGPT to teach me.
Please note that the variable symbols differ from those in the text.

First, calculate the standard error, then calculate the $$${t}$$$ value.

📊 Standard error of the estimator $${\hat{\beta}_j}$$ for the $$${j}$$$-th partial regression coefficient

[Explanation of variables and symbols]

📊 $$${t}$$$ value of the estimator $${\hat{\beta}_j}$$ for the $$${j}$$$-th partial regression coefficient

🔢 Calculation of the test statistic $$${t}$$$ value using the formula
Let's apply the formula above to the orientation data and calculate the $$${t}$$$ value step-by-step.

Extract the explanatory and objective variables from the orientation data, and add a constant term to the explanatory variables.

## 設定と準備
# 説明変数
X = data1[VARS].values
# 目的変数
y = data1[TARGET].values
# 標本サイズN、説明変数の数p
N, p = X.shape
# 説明変数の最初の列に定数項を追加
X_const = np.column_stack([np.ones(N), X])

Calculate the partial regression coefficients.
(I'll calculate the partial regression coefficients while I'm at it)

## 偏回帰係数の推定
beta_hat = np.linalg.inv(X_const.T @ X_const) @ X_const.T @ y
beta_hat

[Execution results]
These are the estimated values for the intercept, temperature, and pressure partial regression coefficients.

Calculate the standard error of the partial regression coefficients.
I'll use the formula!

## 偏回帰係数の標準誤差の推定

# 誤差分散の推定値
sigma2_hat = sum(result1_sm.resid**2) / (N - p - 1)
# 説明変数の行列積の逆行列
XX_inv = np.linalg.inv(X_const.T @ X_const)
# 偏回帰係数の標準誤差の推定
se_beta_hat = np.sqrt(np.diag(sigma2_hat * XX_inv))
# 結果の表示
print(se_beta_hat)

[Execution results]
These are the standard errors for the intercept, temperature, and pressure.
They match the regression analysis results from statsmodels!

Calculate the $$${t}$$$ value and $$${p}$$$ value (two-tailed test) for the partial regression coefficients.

## t値、p値の算出

# t値の算出
t_value = beta_hat / se_beta_hat
print('t値: ', t_value)

# t値のp値の算出
p_value = stats.t.sf(abs(t_value), df=N-p-1)*2 # 両側
print('p値: ', p_value)

[Execution results]
Displayed in the order of intercept, temperature, and pressure.
They match the regression analysis results from statsmodels!

Let's summarize each value calculated using the formula in a table, just like statsmodels does.

## 重回帰分析のサマリー表の作成
pd.DataFrame(
    {'係数': beta_hat, '標準誤差': se_beta_hat, 't値': t_value, 'p値': p_value},
    index=['切片'] + VARS
).round(3)

[Execution results]
It looks pretty good when summarized in a table!

Now, let's proceed to the hypothesis judgment.

③ Determining whether the null hypothesis can be rejected
To state the conclusion first...
For partial regression coefficients other than the intercept, the $${p}$$ value $${\leq}$$ significance level of $${5\%}$$.
At a significance level of $${5\%}$$, the null hypothesis is rejected, and it can be said that the estimated values of the partial regression coefficients other than the intercept are significant.

So, "both temperature and pressure affect the degree of orientation."
That's great!

Let's plot the relationship between the significance level and the rejection region from the text, Figure 2.6.2 "In the case of temperature" and Figure 2.6.3 "In the case of pressure".

### 温度と圧力の係数の検定 p.67 図2.6.2, 2.6.3

## 設定と準備
alpha = 0.05                                 # 有意水準（両側検定）
left, right = -4, 4                          # グラフのx軸の両端
t_vals = result1_sm.tvalues.iloc[1:].values  # t値の取得
p_vals = result1_sm.pvalues.iloc[1:].values  # p値の取得
res = result1_sm                             # statsmodelsの回帰分析結果

## 描画
# 描画領域の設定
fig, axes = plt.subplots(1, 2, figsize=(10, 4), tight_layout=True)
# 係数ごとに描画を繰り返し処理
for t_val, p_val, col, ax in zip(t_vals, p_vals, VARS, axes.flat):
    ## 設定と準備
    # 基本の色
    color = 'tab:blue'
    # 標本サイズn、説明変数の数pの取得
    n, p = int(res.nobs), int(res.df_model)
    # t分布の設定
    t_dist = stats.t(df=n-p-1)
    # 有意水準alphaの両側検定の棄却限界値
    c_val = t_dist.isf(q=alpha/2)
    # t分布の確率密度関数を算出するためのx軸の値
    x_val1 = np.linspace(left, right)   # 確率密度関数用
    x_val2 = np.linspace(left, -c_val)  # 下側棄却域の塗りつぶし用
    x_val3 = np.linspace(c_val, right)  # 上側棄却域の塗りつぶし用

    ## 描画
    # t分布の確率密度関数の描画
    ax.plot(x_val1, t_dist.pdf(x_val1), color=color, label='$t$分布')
    # 下側棄却域の塗りつぶし描画
    ax.fill_between(x_val2, t_dist.pdf(x_val2), color=color, alpha=0.2,
                    label='棄却域')
    # 上側棄却域の塗りつぶし描画
    ax.fill_between(x_val3, t_dist.pdf(x_val3), color=color, alpha=0.2)
    # t値の垂直点線（赤）の描画
    ax.axvline(t_val, color='tab:red', ls='--', label='$t$値')
    # 棄却限界値の垂直点線（黒）の描画
    ax.axvline(-c_val, color='grey', ls='--', label='棄却限界値')
    ax.axvline(c_val, color='grey', ls='--')
    # 修飾
    ax.set_title(f'【{col}】の両側検定の可視化\n'
                 f'自由度{n-p-1}の$t$分布'''
                 f'\n棄却限界値:{c_val:.4f}, $t$ 値:{t_val:.4f} \n'
                 f'有意水準{alpha:.1%}, $p$ 値={p_val:.4f}')
    ax.set_xticks([-c_val, 0, c_val])
    ax.legend(loc='upper left')
plt.show()

[Execution Result]

The blue curve is the probability density function of the $${t}$$ distribution with $${7}$$ degrees of freedom.
The blue regions at both ends of the curve are the rejection regions for a significance level of $${5\%}$$.
The probability (area) of the blue region is $${0.05}$$.
The red vertical dotted line is the $${t}$$ value.
The $${t}$$ values for both explanatory variables are located in the blue region, so the null hypothesis is rejected.

🛸 A little detour: Charts recommended by ChatGPT 🛸

I wanted to know more examples of visualizing multiple regression analysis, so I asked ChatGPT.

The answer I got was...

Let's work on number 4!
Numbers 1 and 2 will be introduced in the next article.
Number 3 was introduced in the previous article.

📈 Volcano Plot
Let's take a quick look at the overview and how to read a volcano plot based on ChatGPT's explanation!

🌋 Principle of Volcano Plot

1. X-axis: Effect size (partial regression coefficient or standardized partial regression coefficient)

   • Indicates how much a variable affects the objective variable by the "magnitude of the value."

2. Y-axis:$${−\log⁡10}$$($${p}$$ value)

   • Represents the smallness of the $${p}$$ value (level of significance) using a logarithm with base 10.

   • For example, $${p\ value=0.01 \rightarrow −\log⁡10(0.01)=2}$$, $${p\ value=0.001 \rightarrow 3}$$, so you can intuitively understand the height that "this point is very significant."

🔍️ How to read

• Points that are far to the right and high up
  → Variables with a large positive effect and are also statistically significant

• Points that are far to the left and high up
  → Variables with a large negative effect and are significant

• Points that are not far to either the left or right, or are located low down
  → Variables where the effect size is small or it is highly likely that they do not meet the significance level ($${p}$$ value is large)

💡 Key Points

• Horizontal spread conveys the "magnitude of the effect size," while vertical height simultaneously conveys the "smallness of the p-value (significance)," allowing you to compare both the "influence" and "reliability" of each variable at a glance.

• The name comes from the idea that if you liken it to a volcano, the variables that pop out on both sides of the ⛰️ look like peaks with "strong significant effects"😊

Using this chart, you can intuitively grasp which explanatory variables are having a "large" and "certain" impact!

We will draw a volcano plot with the horizontal axis (effect) as the standardized partial regression coefficient and the vertical axis (significance) as $${−\log⁡10}$$($${p}$$ value).

### Volcano プロット

## 設定と準備
# 偏回帰係数
betas = result1_sm.params[1:]
# 説明変数の標準偏差、p値を取得
std_X = data1[VARS].std(ddof=1)
# 目的変数の標準偏差
std_y = data1[TARGET].std(ddof=1)
# 偏回帰係数のp値
pvals = result1_sm.pvalues[1:]

## 計算
# 標準化偏回帰係数を計算
beta_star = betas * std_X / std_y
# -log10(p値)の計算
neglogp = -np.log10(pvals)

## 描画
# 描画領域の設定
plt.figure(figsize=(6, 4))
# 標準化偏回帰係数と-log10(p値)の散布図
plt.scatter(beta_star, neglogp, s=70, color='tab:blue')
# 各点に変数名を表示
for var, x_val, y_val in zip(beta_star.index, beta_star.values, neglogp.values):
    plt.text(x_val, y_val, var, fontsize=12, ha='left', va='bottom')
# 有意水準0.05の水平線の描画
plt.axhline(-np.log10(0.05), color='tab:red', ls='--', label='有意水準 5%')
# 修飾
plt.xlim(0.4, 0.6)
plt.ylim(1, 2)
plt.xlabel('標準化偏回帰係数 ($\\beta^*$)', fontsize=12)
plt.ylabel('-log10($p$値)', fontsize=12)
plt.title('Volcanoプロット (標準化偏回帰係数 vs -log10($p$値))')
plt.grid(True)
plt.legend(loc='upper left')
plt.tight_layout()
plt.show()

【Execution Results】

It seems that both the explanatory variables "temperature" and "pressure" have large positive effects and are statistically significant variables.

Comparing the two variables, temperature is

• 【Influence】Large standardized partial regression coefficient
  ⇒ Large influence on the objective variable

• 【Significance】Large $${−\log⁡10}$$($${p}$$ value)
  ⇒ Small $${p}$$ value

it seems we can say.

ChatGPT will conclude the article.
This time, using the analogy of daily growth.

📘 A word from ChatGPT:

That is all for this session of copying code.

Series Articles

Next Article

Previous Article

Table of Contents

Blog Introduction

I am writing a series of seven articles on note.
Please take a look!

1. Relaxed Statistics

This is a blog that roughly explores probability and statistics using the Statistical Test Grade 2 problem collection as a guide.
Feel free to read it as if we were just chatting. Please take a look.
It corresponds to the CBT-compatible version of the official Statistical Test Grade 2 problem collection.
Sample code for Python and EXCEL is also available.

2. Experiment! Fun Bayesian Modeling 1 & 2 with PyMC Ver. 5

I will draw and analyze the Bayesian models used in psychological research from the books "Fun Bayesian Modeling" and "Fun Bayesian Modeling 2" using PyMC Ver. 5.
Like this book, many Bayesian models are written in R language + Stan.
I will strive to explore the possibilities of PyMC and make Bayesian modeling easy to practice.
Since these are familiar and easy-to-visualize themes, please try running them with PyMC and let's enjoy it together!

3. Experiment! Bayesian Modeling from Iwanami Data Science 1 with PyMC Ver. 5

I will draw and analyze the Bayesian models by four Bayesians from the book "Experiment! Iwanami Data Science Vol. 1" using PyMC Ver. 5.
This book is a great resource for learning the basics of Bayesian programming.
I feel like I've become friends with Bayesian methods by happily running PyMC models.
Everyone, please try running them with PyMC, and let's play and learn together!

4. Fun Copying: Bayesian, Python, etc.

I will blog about the results of my "book copying activities" for Bayesian, Python, and others.
I am mainly working on translating them into Python.
I hope this will serve as sample code for fellow learners who are also copying code 🍀

5. Introduction to Time Series Analysis for Psychology with R and Stan, using Python and PyMC Ver. 5

I will practice the time series analysis from the book "Introduction to Time Series Analysis for Psychology with R and Stan" using Python and PyMC Ver. 5.
This book is packed with themes on time series analysis!
I realized the depth of time series analysis.
I will enjoy learning time series analysis with my favorite language, Python.

6. Writing about Data Science-like things

I write columns on statistics, data analysis, AI, machine learning, and Python on an irregular basis.
There are many articles related to statistics and data science books.
Series such as "Statistics," "Python," "Mathematics and Python," and "R" have been created.

7. Python Machine Learning Programming Practice Journal

I wrote articles about my various thoughts when studying the book "Python Machine Learning Programming: PyTorch & scikit-learn Edition."
This book is a textbook for scikit-learn and PyTorch.
Please feel free to try it out if you like.

Thank you very much for reading until the end.