9

Python では、 を使用して r と関連する p 値を計算する方法を知ってscipy.stats.pearsonrいますが、r の信頼区間を計算する方法を見つけることができません。これはどのように行われますか?助けてくれてありがとう:)

4

5 に答える 5

8

[1] によると、Pearson r を直接使用した信頼区間の計算は、正規分布していないため複雑です。次の手順が必要です。

  1. r を z' に変換し、
  2. z' 信頼区間を計算します。z' のサンプリング分布はほぼ正規分布しており、標準誤差は 1/sqrt(n-3) です。
  3. 信頼区間を r に変換します。

サンプルコードを次に示します。

def r_to_z(r):
    return math.log((1 + r) / (1 - r)) / 2.0

def z_to_r(z):
    e = math.exp(2 * z)
    return((e - 1) / (e + 1))

def r_confidence_interval(r, alpha, n):
    z = r_to_z(r)
    se = 1.0 / math.sqrt(n - 3)
    z_crit = stats.norm.ppf(1 - alpha/2)  # 2-tailed z critical value

    lo = z - z_crit * se
    hi = z + z_crit * se

    # Return a sequence
    return (z_to_r(lo), z_to_r(hi))

参照:

  1. http://onlinestatbook.com/2/estimation/correlation_ci.html
于 2017-09-08T03:35:14.427 に答える
1

これは、フィッシャー変換(二変量正規性などを前提とする)ではなく、ブートストラップを使用して信頼区間を計算するソリューションで、この回答から借りています。

import numpy as np


def pearsonr_ci(x, y, ci=95, n_boots=10000):
    x = np.asarray(x)
    y = np.asarray(y)
    
   # (n_boots, n_observations) paired arrays
    rand_ixs = np.random.randint(0, x.shape[0], size=(n_boots, x.shape[0]))
    x_boots = x[rand_ixs]
    y_boots = y[rand_ixs]
    
    # differences from mean
    x_mdiffs = x_boots - x_boots.mean(axis=1)[:, None]
    y_mdiffs = y_boots - y_boots.mean(axis=1)[:, None]
    
    # sums of squares
    x_ss = np.einsum('ij, ij -> i', x_mdiffs, x_mdiffs)
    y_ss = np.einsum('ij, ij -> i', y_mdiffs, y_mdiffs)
    
    # pearson correlations
    r_boots = np.einsum('ij, ij -> i', x_mdiffs, y_mdiffs) / np.sqrt(x_ss * y_ss)
    
    # upper and lower bounds for confidence interval
    ci_low = np.percentile(r_boots, (100 - ci) / 2)
    ci_high = np.percentile(r_boots, (ci + 100) / 2)
    return ci_low, ci_high
于 2020-10-14T21:55:19.947 に答える