SSE-DP-2024-3

t 統計量の分布が双峰型となる場合

国友直人・西颯人・薛玉傑

統計数理研究所

2024 年 7 月
2024 年 9 月(改訂版)

SSE-DP(ディスカッションペーパー・シリーズ)は以下のサイトから無料で入手可能です。
https://stat-expert.ism.ac.jp/training/discussionpaper/
このディスカッション・ペーパーは、関係者の討論に資するための未定稿の段階にある草稿である。著者の承諾なしに引用・複写することは差し控えられたい。

SSE-DP-2024-3
On Bimodality of t-StatisticbyNaoto Kunitomo, Hayato Nishi and Xue Yujie

The Institute of Statistical Mathematics

July 2024
September 2024 (Revised)（Summary）
Many statisticians use the t-distribution for t-statistic under the assumptionof normal population. The distribution of t-statistics can be bimodal when
the population distribution is Cauchy. Although there were some studies
on this issue by Fioro, Hajivalliliou and Phillps (2010, Econometric Journal), andLogan, Mallows, Rice and Shepp (1973, Annals of Probability), theiranalysis are difficult to be understood for applied statisticians. We discuss
the mathematical foundations for applications and also report somesimulation results when the population is stable distributions.

統計エキスパート演習 (統計基礎を巡る話題 2024):「t 統計量の分布が双峰型となる場合」 12024 年 7 月 26 日
2024 年 9 月 18 日 (改訂) 2
国友直人3 ・西颯人4 ・薛玉傑5

要約
統計数理研究所が推進している統計エキスパート養成プログラムの意見交換サイト（slack）
上で「t-統計量の標本分布が双峰型になるシミュレーション」と云う（三輪哲久特任教授が掲示した）R プログラムを巡って、「t 統計量の標本分布」について活発な議論が行われた。コーシー分布を含む安定分布などについてのシミュレーションの結果、数理的基礎および関連する統計的問題について関心を抱いた参加者の考察を報告する。また議論の内容の理解に資すると思われる事項を注として述べるとともに付論 A として利用した計算プ
ログラム、付論 B, 付論 C としてこれまであまり複素関数論・確率論などを学ぶ機会がな
かった方々のために図を掲載、基礎的事項を引用、応用統計家にとって有用と考えられる話題を解説した。

はじめに（国友）

1

統計検定 2 級の教科書「統計学基礎」3 章では多くの教科書と同様に t 統計量についての標準的な説明がある。正規分布にしたがう母集団からの標本平均は母集団の分散が未知の時には標本から計算される分散で基準化すると、その分布は正規分布ではなく t 分布にしたがうことが説明されている。また統計学の歴史を紐解くとダブリン (アイルランド) にあるビール会社ギネス社の技術職員だった William Gosset による t 統計量の分布の導出(Biometrika, 1905) がその後の近代的な統計学の展開の鍵として (こうした業績の評価はR.A. Fisher から始まると思われるが)Wikipedia を含め、説明されている。

例えば日本の統計家の間で良く引用されている竹村彰通「現代数理統計学」では互いに独立で同一の分布に従う確率変数列 X1 , X1 , · · · , Xn について、もし母集団が正規分布∑
N (µ, σ 2 ) であるなら標本平均 X̄n = (1/n) ni=1 Xi について µ を既知としたときの Tn 統

1

統計エキスパート slack 上で行われた統計的推測の基礎についての素朴な疑問や議論の中から統計エキ
スパート養成事業にとり有用と思われる内容をまとめたメモである。三輪哲久特任教授の問題提起 (R プロ
グラムと図 1) より議論が始まったので同教授の許可のもとに原プログラムも掲示しておく。
2
SSE-DP-2024-3(統計エキスパート DP) の改訂版
https://stat-expert.ism.ac.jp/wp/wp-content/uploads/2024/07/SSE-DP-2024-3.pdf
3
統計数理研究所
4
東京大学工学系研究科都市工学専攻
5
統計数理研究所1

計量√n[X̄n − µ]Tn = √ ∑n2i=1 (Xi −X̄n )

(1.1)

n−1

が自由度 n − 1 の t 分布にしたがうことが説明されている。ここで自由度 m の t 分布の密度関数は基準化定数 c(m) とおくと f (x) = cm [1 + x2 /m]−(m+1)/2 (−∞ < x < ∞) で与えられる。続いて同書 page 81 では正規母集団でなくても標本不偏分散は∑n2i=1 (Xi − X̄n )
(1.2)
Vn =
n−1
は n が十分に大きいとp

Vn −→ σ 2

(1.3)

と確率収束 (convergence in probability) する。したがって、母集団か必ずしも正規分布でなくても標本平均について中心極限定理 (Central Limit Theorem, CLT) より√wn[X̄n − µ] −→ N (0, σ 2 )(1.4)と分布の意味での収束、分布収束あるいは弱収束 (weak convergence) するので、統計量
Tn についても標準正規分布へ収束してw

Tn −→ N (0, 1)(1.5)となることも説明されている。こうした日本の統計家の間で標準的とされている書籍における説明は、（CLT が成り立つ条件などを巡る古典的確率論について内容の正確な理解はさておき）今日でも統計的方法を利用する応用家の一つのよりどころとなっていることは間違いない。

こうした統計量 Tn の分布についての説明があるにも拘わらず、母集団が正規分布でない場合に t 統計量の分布が歪むこと、例えば二山型になることはあり得ないのだろうか？
こうした素朴ではあるが統計分析での基本的な問題について統計エキスパート養成事業のメンターの一人、三輪哲久教授は短い R プログラムと次のような事実を同事業の slack
上で指摘した。
「母集団がコーシー分布 Cauchy(0,1) にしたがう場合、シミュレーションにより Tn 統
計量の分布は図１の様に双峰型になる！」6
近年ではフリーの統計計算ソフトウエア R、Python などによりこの種のシミュレーショ
ンを実際に行うことはかなり容易である。この指摘を受けて様々な実験が可能であるが、分析対象を標準的な i.i.d. 確率変数としても、次のような幾つかの基本的な疑問が生じる。

(i) Cauchy 分布は自由度 1 の t 分布である。母集団分布が仮に自由度 n の t 分布である場
合、Tn 統計量の分布はどうなるだろうか？自由度は必ずしも正整数でない実数でも定義は可能である。
6

三輪教授の R プログラムは強化の上で付論 A に掲載しておいた。

2

(ii) Cauchy 分布は安定分布族 (stable distribution) の母数 α = 1 の場合であることが知られている。（定義は付論を参照のこと。）安定分布族は 0 < α ≤ 2 のパラメターを持ち、
α = 2 が正規分布となることが知られている。母集団分布が仮に安定分布である場合、Tn
統計量の分布はどうなるだろうか？なお、安定分布族については付論 B でその定義と性質を簡単に述べておくが、統計学で用いられている安定分布にしたがう標準的な乱数発生法については４節を参照されたい。
(iii) 統計的データ分析における応用上では Tn 統計量の分布が t 分布や正規分布とはかな
り異なる場合、すなわち漸近正規性が成り立たない場合はあり得るのだろうか？その場合にはどの様な統計的方法が有用なのだろうか？

2

t 統計量の分布が双峰型となる実験 (西)三輪教授が作成した R によるシミュレーション・プログラムを修正、Python プログラム
による幾つかのシミュレーション実験を行った。念のために作成したプログラムは付論 Aに掲載して添付しておいた。

2.1. (i) についてのシミュレーション実験：
t 分布の自由度とサンプルサイズをいろいろ変えて（n = 10, 30, 100, 1000）シミュレーショ
ンを実行してみると、母集団が df=1（コーシー分布）だけは双峰になり、サンプルサイ
ズを大きくしても双峰性は維持される。すなわちサンプルサイズの影響は小さい。平均値への影響が大きい外れ値の発生確率を考えると、上側と下側で均等に発生する確率よりどちらかに偏る確率の方が高いからこうなるのかもしれないとも考えられる。ここで「平均値への影響が大きい外れ値」をどう定義すればよいのか検討の余地がある。
2.2. (ii) についてのシミュレーション実験
安定分布（stable distribution）についてもサンプルサイズ（n = 10, 30, 100, 1000）として母集団が安定分布となる場合のシミュレーション実験を行ってみた。なお一般の安定分
布の場合には分布の歪みも考えられるが、今回は対称分布に限定して実験を行った。母数
α = 1.5 以上だと谷が消え、単峰性があるようである。ただし安定分布の非対称性を示す指標である β = 0 に固定、ここでは 1 ≤ α ≤ 1.5 の結果を示しておいた。この実験における結果 α = 1.5 についてはこうした現象について言及している文献は今のところ見つかっていない。更に検討する価値はあるだろう。

3

100005000
0Frequency

15000

Histogram of t.cauchy

-6

-4

-2

0

2

4

6

t.cauchy

図 1: t 統計量の分布 (シミュレーション)

4図 2: t 統計量の分布 (シミュレーション)

5図 3: t 統計量の分布 (シミュレーション)

6

漸近分布の数理について（薛）3i.i.d. 確率変数列が Cauchy 分布にしたがう場合の t 統計量の分布については既に [3] が議論している。この研究はファイナンス分野における裾の厚い分布の統計分析における応用問題の議論から問題を考察したようである。ただし [3] が利用している数理的議論はそのまま直ぐに理解できるとは言い難いので、特殊関数についてよく知られている文献 [2] および確率論の文献 [4] に基づきなるべく分かりやすい解説を試みる。
2 階常微分方程式 (second-order ordinary diﬀerential equation)
xy ′′ + (c − x)y ′ − ay = 0(3.1)

y = b1 Φ(a, c; x) + b2 Ψ(a, c; x),

(3.2)

の解の形式
を考察する。ここで Φ(a, c; x) = 1 F1 (a; c; x) は第一種の合流型超幾何関数 (confluent hy∫ ∞ −xt a−11pergeometric function of first kind) と呼ばれているが、Ψ(a, c; x) := Γ(a)e t (1 +0
c−a−1t)
dt for Re a, Re x > 0 は第二種の合流型超幾何関数 (confluent hypergeometric function of second kind) と呼ばれている。これらの関数には

Ψ(a, c; x) =Γ(c − 1) 1−cΓ(1 − c)Φ(a, c; x) +x Φ(a − c + 1, 2 − c; x),Γ(a − c + 1)
Γ(a)

(3.3)

という (3.3) 関係がある7 。なお一般の超幾何関数 (generalized hypergeometric function) と∑∞ (a ) ···(a ) k
は p Fq (a1 ...., ap ; b1 , ..., bq ; z) = k=0 (b11 )kk ···(bqp)kk zk! により定める。ここで記号 (a)k = a(a +1) · · · (a + k − 1), a0 = 1 により定める。
超幾何関数を利用して Cauchy 分布からの Tn -統計量の性質を調べよう。X1 , . . . , XT を
1
Cauchy (0, 1) からの i.i.d. 標本とする。密度関数は f (x) = π(1+x2 ) (−∞ < x < ∞) である
1
2
2
が、limx→∞ x f (x) = limx→∞ x f (−x) = π を満足している。
∑T
∑T
Xk
2
2
2
2
ここで UT = k=1,V=T
i=1 Xi /T および ST = UT /VT として、統計量 VT の特T性関数を考察する。i.i.d. にしたがう確率変数 X とすると

cfVT2 (t) = {cfX 2 /T 2 (t)}T= { E exp(i|X|2 t/T 2 )}T∫ ∞exp(itx2 /T 2 ) dx T={}
π(1 + x2 )−∞∫ ∞exp(itr/T 2 ) dr T={}
πr1/2 (1 + r)
0
= {A(t)}T
7

(say)

(3.4)

4 節の注を参照されたい。特殊関数については文献 [2] が詳しい。(3.3) 式は p.257 の (7) 式に与えられている。ただし証明にはガンマ関数に基づく二つの関数 Ψ(·) と Φ(·) の複素積分表示から注３で述べるような（複素数上での）Γ 関数の極と留数 (residue) の定理を利用する必要がある。特殊関数についての正確な評価が必要となるのでここでは省略する。

7

と表現できる。ここで (3.3) より

A(t) = lim+
ε→0

Γ( 12 ) 1 1
Ψ( , ; −it/T 2 + ε),
π
2 2

(3.5)

であるから

Γ( 21 ) 1
Γ(− 21 ) (−it)1/2
1 1
3−2 (−it)1/2
1
2
2
A(t) =[Γ( )Φ( , ; −it/T )+
Φ(1,
;
−it/T
)]
=
1+
+o( ).
1
1/2
π
2
2 2
T
2
π
T
T
Γ( 2 )(3.6)∫ ∞ z−1∑∞n1となる。ここで Γ(z) = 1 t et dt + n=0 (−1)
は拡張されたガンマ関数 Γ であり、通n! n+z常のガンマ関数 Γ(z) の定義域 Re z > 0 から複素平面上 (ゼロおよび負の整数を除く) に連続的（解析的）に拡張した関数を意味する8 。
1/2
さらに極限をとると limT →∞ cfVT2 (t) = exp[ π−2] となる.
1/2 (−it)
ここで −it を極座標表示すると −it = |t| exp(−isgn(t)π/2) であるから(−it)1/2 = |t|1/2 exp(−isgn(t)π/4) = |t|1/2 cos(π/4)(1 − isgn(t) tan(π/4)).

(3.7)

したがって Lévy’s 連続性定理 (continuity theorem) より統計量 VT2 の極限分布の特性関数
1/2
は exp{ π−2cos(π/4)(1 − isgn(t) tan(π/4))]} で与えられる。この特性関数の形から
1/2 [|t|α = 1/2 の安定分布であることが分かる。
次に (UT , VT2 ) の特性関数を考察する。まずcfUT ,VT2 (s, t) = {cfX/T,X 2 /T 2 (s, t)}T .
∫ ∞

2

(3.8)

2

eisx/T +itx /T
cfX/T,X 2 /T 2 (s, t) =dx.
2−∞ π(1 + x )(3.9)となることを利用する。ここで変換 x → (r, h) (ただし r = x2 , h = sgn(x)）を利用する.
このとき
∫ ∞
∫
√
exp(itr2 /T 2 )1
cfX/T,X 2 /T 2 (s, t) =exp(is rh/T ) dh dr.
(3.10)2π(1 + x )0h 2となる。ここで Taylor 展開を ez に施すと次の結果が得られる。

∫√(is/T ) rh

1/2

e

dh =

h

∞
∑
((is/T )2 /4)k rk
k=0

k!( 12 )k

.

したがって

∫ ∞
cfX/T,X 2 /T 2 (s, t) =0∞
exp(itr2 /T 2 ) ∑ ((is/T )2 /4)k rkdrπ(1 + x2 ) k=0
k!( 21 )k1
1
1 ∑ ((is/T )2 /4)k Γ(k + 21 )
Ψ(k + , k + , −it/T 2 )
=
1
π k=0
2
2
k!( 2 )k∞8

(3.11)

ここでは Erdélyi [2], Chapter 1 の表現を採用した。正実数に対するガンマ関数を複素関数へ拡張する
ときには幾つかの等価な定義が存在する。複素解析についての標準的文献としては Ahlfors, L.V. (1979),
Complex Analysis, 3rd Edition, McGraw-Hill, Inc., 高橋礼司 (1989), 「複素解析」, 東京大学出版会、読みやすい書籍としては志賀浩二 (1989), 「複素数 30 講」, 朝倉書店 などがある。

8

となる。(3.3) を用いて Ψ(k + 12 , k + 12 , −it/T 2 ) に Γ( 12 − k)Φ(k + 21 , k + 12 ; −it/T 2 ) +Γ(k− 21 )π1
12(−it/T 2 )1/2−k Φ(1, 32 − k; −it/T 2 ), Γ( 12 − k) = (−1)k Γ(k+1 , Φ(k + 2 , k + 2 ; −it/T ) =Γ( 12 +k))2
−it/T 2eを代入するとcfX/T,X 2 /T 2 (s, t) =∞∑(−(is/T )2 /4)kk=0∞∑

k!( 21 )ke−it/T2((is/T )2 /4)k Γ(k − 12 )1
3+(−it/T 2 )1/2−k Φ(1, − k; −it/T 2 ) (3.12)
1
π k=0
2
k!( 2 )kとなる。ここで Γ(k − 12 ) = Γ(− 21 )(− 12 )k となることに注意すると,∞
Γ(− 12 ) ∑ ((is/T )2 /4)k (− 21 )kcfX/T,X 2 /T 2 (s, t) = 1+o(1/T )+(−it/T 2 )1/2−k (1+o(1/T ))π k=0
k!( 21 )k∞Γ(− 21 ) ∑ ((s2 /4it)/4)k (− 12 )k=1+()(−it/T 2 )1/2 + o(1/T )
1
π
k!(
)
k
2
k=0

=1+

Γ(− 12 )1 1
Φ(− , ; s2 /4it)(−it/T 2 )1/2 + o(1/T ) (3.13)π2 2が得られる。したがって極限では−2
1 1
lim cfUT ,VT2 (s, t) = exp[ 1/2 Φ(− , ; s2 /4it)(−it)1/2 ].
T →∞
π
2 2

(3.14)となるので、二次元確率変数 (UT , VT2 ) の極限分布の特性関数が得られた。そこで次にS = U/V の極限分布を導こう。
再び特性関数−2
1 1
E exp(iU s + iV 2 t) = exp[ 1/2 Φ(− , ; s2 /4it)(−it)1/2 ] = exp[B(s, t)] (say).
π
2 2(3.15)

を考察する。ここで正実数 a,b をとり t = is2 a および t = is2 b とおき、0 < s < ∞ で積分
すると
∫ ∞
∫ ∞
1 iU s −V 2 s2 a
1 sB(1,ia)−V 2 s2 b
E e {e−e} ds =(e− esB(1,ib) ) ds
(3.16)ss00となる。さらに V s = t とおくと∫ ∞∫ ∞1 i(U/V )t −t2 a
1 sB(1,ia)
−t2 b
E e
{e
− e } dt =(e
− esB(1,ib) ) ds
t
s
0
0
∫ ∞
∫ ∞
1
1 sB(1,ia)
−t2 a
−t2 b
ϕ(t){e− e } dt =(e
− esB(1,ib) ) ds,ts00となる。ここで ϕ(t) は S の特性関数である。(3.17) の両辺を a で微分すると、∫ ∞∫ ∞B ′ (1, ia)−t2 a′tϕ(t)edt = −B (1, ia)esB(1,ia) ds =.
B(1, ia)00
さらに (πa)−1/2 e−s /4a を乗じて、a について 0-∞ で積分し、恒等式∫ ∞texp[−s2 /4a − at2 ]a−1/2 da = exp(−st),s > 0, t > 0,1/2
π0(3.17)(3.18)

2

9(3.19)

を利用すると s > 0 に対して
∫ ∞∫ ∞B ′ (1, ia)2−stϕ(t)e dt =(πa)−1/2 e−s /4ada.
B(1, ia)00(3.20)

が得られる。
確率変数 S の密度関数を g(y) とする. すなわち、特性関数の逆変換 (実部)∫1 ∞ϕ(t)e−st dt
g(y) = lim Re
s→iy
π 0

(3.21)∫∞′ (1,ia)2
da を s = iy としで与えられる。ここで [3] の (4.35) を利用して 0 (πa)−1/2 e−s /4a BB(1,ia)て連続的に変形すると、積分値は∫ ∞1g(y) = lim Re
s→iy
π

−st

ϕ(t)e

∫ ∞/i 12dt = Re1ey 2 t2 /2(3.22)

00と表現できる9 。ただし D̃(t) = ( π23 ) 2 {1 + tet /2iπ
さらに変換 t → re− 4 を用いると∫ ∞/i 122 2

ey t /2 D̃(t) dt

D̃(t) dt =

0

2

∫ ∞

∫ t −u2 /2edu}−1 で与えられる。
0e−y r i/2 D̃(re− 4 )e− 4 dr
iπ

2 2

iπ

(3.23)

0と表現できる。
1
ここで [4] の (iii) より （y → 1 のとき）, g(y) ∼ π12 log |1−y|となる。また Logan et al. [4]の補題 A , 補題 B を利用すると y → ∞ のとき g(y) ∼ ( π2 )1/2 τ0 e−1/2y τ0 となる。ただしD̃(−iτ0 ) = ∞ である.
2 2

最後に t 統計量の漸近分布を考察しよう。t 統計量 tT =UT

1

(VT2 − T1 UT2 ) 2とする. UT = X̄ ∼Cauchy(0, 1) のとき VT2 − UT2 /T = VT2 − Op (1/T ) である. したがって tT = ST + Op (1/T )となる. すなわち t 統計量の極限分布は ST の極限である S. の分布に一致することが分かる。
注 1 Logan et al. [4] の (4.35) について :Logan et al. [4] より B(1, ia) は µ = 0 としたときの ψ̃(a) （Logan et al. [4] の (4.41))、す
なわち
∫ ∞
1
2
[e−x a eix − 1] 2 dx
B(1, ia) =(3.24)x−∞である。したがって′∫ ∞B (1, ia) =−∞−e−x a eix dx.
2(3.25)となる。ここでa = 1/(2b2 ) , x = ξb
9(3.26)

この変形は Cauchy の積分定理による。Logan et al. [4] に極の位置に関する正当化が説明されているがかなり技術的議論なのでここでは省略する。なお（3.22）の右辺は複素変数 t を極座標 (r, θ) で表示すると、
i1/2 = [exp(iπ/2)]1/2 = exp (iπ/4) より t = r ∗ exp(−iπ/4) での実数 r の範囲は (0, ∞)、偏角 θ = −π/4 に設定という意味である。ここで dt ＝ exp(−iπ/4)dr とも表現できる。

10

を代入すると右辺 (3.20) は∫ ∞′−1/2 −s2 /4a B (1, ia)

(πa)

e

B(1, ia)

0

∫ ∞√
da =
0

∫ ∞

=π

∫∞
2
−e−ξ /2 eiξb dξ
2 −s2 b2 /2
−∞
∫∞
e
db
π
[e−ξ2 /2 eiξb − 1] ξ12 dξ−∞e−s b /2 D̃(b) db,
2 2(3.27)

0となるが、ただし D̃(b) =√∫ ∞ −ξ2 /2 iξb
√
∫∞
e dξ
2
2 ∫
1
−∞ e
を乗じた分子 −∞ e−ξ /2 eiξb dξ
である。
π 3 ∞ [1−e−ξ2 /2 eiξb ] 12 dξ
2π−∞は正規分布の特性関数∫ ∞ξ

e−ξ /2 eiξb dξ =2√22πe−b /2

(3.28)−∞

である。分母は∫ ∞∫ ∞12−ξ 2 /2 iξb 1[1 − ee ] 2 dξ =[1 − e−ξ /2 cos(ξb)] 2 dξ
ξ
ξ
−∞
−∞

∫ ∞ −ξ2 /2
∫ ∞
1∞
be
sin(ξb)2= −[1 − ecos(ξb)] |−∞ +dξ +e−ξ /2 cos(ξb) dξ
ξ
ξ
−∞
−∞
∫ ∞
∫ b
√
2
2
=0+b
e−ξ /2 (cos ξu du) dξ + 2πe−b /2−ξ 2 /2

−∞

0

∫ b√
√
2
2
2πe−u /2 du + 2πe−b /2=b(3.29)

01である。したがって D̃(t) = ( π23 ) 2 {1 + tet /22∫ t −u2 /2edu}−1 が得られる。
0

注 2 : 複素積分による評価について (極と留数の定理など)
2 2
Logan et al. [4] 補題 A より負の虚軸 t 上で関数 ey t /2 D̃(t) は D̃(t) = ∞ となるただ一つの極 (pole) t0 が存在する。そこで閉曲線 C(ρ)、曲線 C1 (ρ) (0 から ρ exp −iπ/2 および第４象限の t0 の周りを修正), C2 (ρ) (ρ exp −iπ/2 から ρ exp −iπ/4)、C3 (ρ exp −iπ/4 から
0) をとる。このとき
∫
∫
∫
∫
2 2
y 2 t2 /2y 2 t2 /2y 2 t2 /2eD̃(t) dt =eD̃(t) dt +eD̃(t) dt +ey t /2 D̃(t) dt
C(ρ)

C1 (ρ)

C2 (ρ)

C3 (ρ)(3.30)∫

2 2となる。Logan et al. [4] 補題 B から ρ∗n → ∞ を C2 (ρ∗ ) ey t /2 D̃(t) dt ((3.30) の右辺) がゼnロに収束して、C2 (ρ∗n ) と C2 (ρ∗n+1 ) の間に唯一の極があるようにとれる。この時の極と経路 C(ρ) を図 4 で表す。
ここで留数の定理 (The residue theorem) を用いると∫exp(y t /2) D̃(t) dt = 2πi
2 2C(ρ∗n )n∑Res(exp(y 2 t2 /2) D̃(t), ti )(3.31)

i=1となるがここで ti は C2 (ρ∗i ) と C2 (ρ∗i−1 ) の間にある極である。したがって n → ∞ のとき(3.30) より∫∫∫
2 2
y 2 t2 /2y 2 t2 /2
ey t /2 D̃(t) dt. (3.32)
eD̃(t) dt + limeD̃(t) dt = lim
limn→∞

C(ρ∗n )n→∞n→∞

C1 (ρ∗n )11

C3 (ρ∗n )図 4: 極 (poles) と C(ρ)
となる10 。さらに (3.31), (3.32), (3.22) より C2 (ρ∗n ) 上の積分をゼロとして経路の符号を考
慮すると、
∫
2 2
ey t /2 D̃(t) dt}(3.33)
g(y)=Real{− limn→∞∫

C3 (ρ∗n )y 2 t2 /2=Re{ limn→∞e
C1 (ρ∗n )D̃(t) dt − 2πi

∞
∑

2 2

Res(ey t /2 D̃(t), ti )}.」

i=1が得られる。虚軸上で関数 ey t /2 D̃(t) の実部は 0、さらに limn→∞
2 2∫
C1 (ρ∗n )ey t /2 D̃(t) dt の2 2

実部 −1/2 Re{2πiRes(ey t /2 D̃(t), t0 } は t0 周辺で修正した経路上での積分値の実部に対応√2 2する。Logan et al. [4] の (5.6) より Res(ey t /2 D̃(t), ti } = −2 2

2 y 2 t2i /2eti for i = 0, 1, 2, ...
π3となるので11g(y) = Re{2πi

∞
∑
i=1

√

2 y2 t2i /2eti } + Re{πiπ3√2 y2 t20 /2
e
t0 }
π3

(3.34)となる。Logan et al. [4] 809 項の議論より |y| > 1 のとき g(y) は有界、さらに数値的検討
2 2
により Logan et al. [4] は y → ±∞ のとき g(y) ∼ ( π2 )1/2 τ0 e−1/2y τ0 となると推測している。なお τ0 は τ0 > 0 , D̃(−iτ0 ) = ∞(虚軸上の極) である。また |y| < 1 では Logan et al.
[4](4.36) の議論 (実数 t をとると D̃(t) = O(exp(−t2 /2)) など) より∫ ∞2 2
g(y) =
ey t /2 D̃(t) dt,(3.35)
010

ここで有限個の極が存在するとき可積分であることが示せるので Logan et al. [4] が議論しているように極限をとることが正当化できると思われるので詳しい議論は省略する。
2 2
11
ここで孤立特異点 ti (i = 1, · · · , n) における留数 (residue) は項 ey ti /2 (2/π 3 )1/2 [−ti /(t − ti )] の係数
2 2ey ti /2 (2/π 3 )1/2 (−ti ) を利用した。

12

t → ∞ のとき
2 2

ey t /2 D̃(t) ∼1 (y2 −1)t2 /2 −1et
π2(3.36)と近似できる。したがって∫ ∞1 (y2 −1)t2 /2 −1g(y) ∼et dt + C
2
1−y π
∫ ∞
1
2
∼ 2
e(y −1)s/2 /s ds + C
2π (1−y)2∫ ∞1 (y2 −1)s/2
2∞(y 2 − 1)/2e(y −1)s/2 ln s ds + C(3.37)∼ 2eln s|(1−y)2 −2π(1−y)2∫∞と近似できる。ただし 0 e−s ln s ds = −γ となることに注意すると (ここで γ は Euler の定数) 近似的に
1g(y) ∼ − 2 log|1 − y|(3.38)π(y → 1 のとき) となる。同様に y → −1 のとき

1g(y) ∼ − 2 log|1 + y|π

(3.39)となることが分かる。

4

幾つかの注意点 (国友)

本稿では統計学ではよく知られた Tn 統計量の標本分布についてあらためて 3 名が検討した結果を報告するが、この間に気が付いた統計エキスパートにとり有用と考えられる幾つかの注意点をまとめて述べておく。

(i ) 安定分布にしたがう確率変数の生成
正規分布にしたがう確率変数列を生成する方法としては Box-Muller 法が統計学では良く知られている。一様分布 U (0, 1) にしたがう互いに独立な確率変数 U1 , U2 より変換 X1 =√√µ + σ −2 log(U1 ) cos 2πU2 , X2 = µ + σ −2 log(U1 ) sin 2πU2 , は独立な N (µ, σ 2 ) にしたがう二つの正規確率変数となる。
安定分布にしたがう確率変数列については Chamber-Mallows-Stuck 法 (Chamber et al.
(1976)) が統計学では標準的である。Θ を一様分布 U (−π/2, π/2)、W を λ = 1 の指数分
布に独立にしたがうとする。
1. 0 < α ≤ 2, α ̸= 1 のときZ=

sin αΘ cos(α − 1)Θ (1−α)/α],[(cos Θ)1/α
W

(4.1)

α = 1 のとき Z = tan Θ とする。このとき Z ∼ S(α, 0, 1), 安定分布 S(α, 0, 1) にしたがう。
2. 一般の (歪みのある) 安定分布を含む場合は、−1 ≤ β ≤ 1, θ0 = arctan(β tan(πα/2))/α(α ̸= 1) とおく。0 < α ≤ 2, α ̸= 1 のときZ=

sin α(θ0 + Θ)cos(αθ0 + (α − 1)Θ (1−α)/α
]
,
[
(1/α)
(cos αθ0 cos Θ)
W
13

(4.2)

α = 1 のときZ=π
W cos Θ2 π[( + βΘ) tan Θ − β log( 2 π
)]
π 2
+
βΘ
2

(4.3)

とする。このとき Z ∼ S(α, β, 1), 安定分布 S(α, β, 1) にしたがう。
（証明は Nolan (2020)3
章に与えられている。）(ii) T n 統計量の分布
竹村彰通「現代数理統計学」など日本の統計学関係者の間では定評がある教科書における中心極限定理 (CLT) の説明は十分とは言えない。古典的な確率論では大数の法則や CLTについて詳しく調べられているが、独立な確率変数列についての CLT では Lindeberg 条
件が基本的であり、よく知られた十分条件として確率変数が二次積率 E[X 2 ] < ∞ を持つことが要請される12 。コーシー分布の密度関数は f (x) = (1/π)(1/(1 + x2 ) (−∞ < x < ∞)であるから期待値 E[|X|] は発散する。自由度 n の t 分布は n = 1 の場合はコーシー分布、
n = 2 の場合は期待値は有限であるが二次積率は発散する。安定分布としては α = 2 は正規分布となるが、α = 1 はコーシー分布に対応している。したがって 1 < α < 2 の場合については幾つかの興味深い問題が生じ得る。CLT が成立するとは限らないが、再生性（和の分布が元の分布に属するという性質）があることから古典的な確率論ではかなり研究されている13 。しかし、コーシー分布からの統計量の分布についてはあまり注目されてこなかったようである。
このレポートの３節において解説しているように数理的には実はかなり難しい議論が必要となることもこの問題が十分に理解されていない原因と考えられる (参考となる事項を以下に追加の注として述べる)。第１節で述べたように Tn の漸近正規性が成り立つ
ためには (i) 分子について CLT が成り立つ, (ii) 分母が一定の正値に確率収束する、ことが必要となる。コーシー分布の場合には位置母数 µ = 0 としておくと (i) コーシー分布∑nからの標本平均 X̄n = (1/n) i=1 Xi は再びコーシー分布にしたがう。（付論 B を参照。）∑(ii)Vn = ni=1 Xi2 /n2 は n −→ ∞ のときある確率変数 V (> 0) に収束する。この (ii) を示
すことは自明ではないので正確に示すことが第３節の最初の課題である。分布を明示的に表現できないので特性関数による表現を利用するが、特性関数と分布関数は 1 対 1 に対応
する（Levy の反転公式14 による）ので、特性関数の極限が存在すれば極限分布が存在することが分かる。コーシー分布の場合には V は α = 1/2 の安定分布となることが示される。次に X̄n と Vn の同時分布の極限が存在することが３節の第２の内容である。こうして得られた同時特性関数から統計量 Tn の分布を導くには有限個とは限らない極 (poles) に留数の定理 (residue theorem) を適用する必要があるようである。したがってこの辺からなぜ問題が簡単でないことも理解できる。さらに最終的な目的である Tn 統計量の極限分∑n2
独立な確率変数列 Xm,n (1 ≤ ∑m ≤ n) が E[Xn,m ] = 0, m=1 E[Xn,m] −→ σ 2 (> 0) 及び任意の ϵ > 0n2に対して Lindeberg 条件 limn→∞ m=1 E[Xn,m ||Xn,m | ≥ ϵ] −→ 0 が成り立つとき S = Xn,1 + · · · + Xn,nについて CLT が成立する (例えば Durrett (1991), Page 98 参照)。Lindeberg 条件は CLT の必要十分に近いが独立性は Martingales
や定常確率過程に緩めることもできる。統計学の普通の教科書に書いてある CLT√
は Xn,m = Xm / n (ただし Xm は独立な確率変数列) と云う特殊な場合が多い。
13
CLT や安定分布などを含め古典的確率論においてよく知られている内容を説明している教科書の例としては Durrett, R. (1991),”’Probability : Theory and Example,”Duxbury を挙げておく。
∫T
14
1 次元確率測度 µ(·), 特性関数 ψ(·) とすると任意の実数 a < b に対し limT →∞ (2π)−1 −T [e−iat −e−ibt ]/(it) × ϕ(t)dt = µ(a, b) + (1/2)µ(a, b) となる。例えば Durrett (1991) Page 78 を参照。
12

14

布は存在はするものの、CLT による正規分布とは様相がかなり異なることが分かる。例
えば極限分布の密度関数はある点において発散しうる、などの論点の説明が３節の内容である。
こうした Cauchy 分布を巡る結果は既に Fioro et al. (2010) の Theorem 2.1 として述べられているが、彼らの証明は不完全と思われ（あるいは分かりにくい）ので、本稿ではその簡単な導出を試みている。ところで、問題を検討する過程で実はこうした結果は Loganet al. (1973) により既に安定分布に関する様々な結果の一部として既に得られていたことが判明した。しかしながら、Logan et al. (1973) の論文は (普通の統計家には) かなり高度な複素解析・特殊関数の議論を多用して数理的結果を導くことに加えて、数学関係者を対象とする確率論の専門学術誌の論文であり、普通の応用統計家にはかなり近寄りがたい感がある。そこで３節の後半では Logan et al. (1973) による結果を可能な範囲でより分かりやすく理解できるように彼らの論文には説明がない幾つかの議論と解釈を展開した。
Logan et al. (1973) の図 5 では α = 1.5 の場合について密度関数は発散しないが僅かに双峰となる場合があるようである。また極限分布を通常の反転公式を利用した積分計算により数値計算を試みると、しばしば数値が発散した。これは極限分布の密度関数が場所により発散しうることに符合すると解釈できるかもしれない。
本稿での考察から結論として得られた Tn 統計量の漸近分布が双峰型になり得る、という論点は少なくともこれまで多くの日本の統計家の間で知られていたとは言い切れない。
応用上の可能性を含めて一つの重要な結果と言えるだろう。

(iii) 安定分布と応用
安定分布の株価への応用可能性を 1963 年15 に指摘したのは物理学者の Mandelbrot であり、その後、一時期は E. Fama に代表されるファイナンス分野の研究者が注目した。その後、こうした流行はかなり下火になったとはいえ、再び金融リスク管理の問題においても確率分布の裾 (稀に起きる事象 rare events) の評価に関連して注目されている。株価と云ったファイナンス分野に限らず工学や理学分野などでの統計データ分析における応用可能性については例えば Nolan(2020) が積極的な立場から展望を行っている。
(iv) 幾つかの残された課題
正規分布やコーシー分布を含む確率分布として安定分布族の重要性は理解できたとしても、その先はあり得るだろうか？確率変数列の和の極限分布としては正規分布やコーシー
分布が挙げられるが、確率論では安定分布族を含む、和についての全ての可能な極限分布を表現する無限分解可能分布 (付論を参照されたい) があり、安定分布は一つの特殊な場合であることが知られている。こうした一般化により統計的データ分析における新たな手がかりが得られる可能性は皆無ではないと思われる。
(v) 関連する応用統計の問題
応用上で表れる母集団の積率を巡る統計的問題にも言及しておこう。単純な例として母集団が N (1/µ, 1)(µ ̸= 0) であることが分かっている場合を考えよう。母集団から i.i.d. 標本
Xi (i = 1, · · · , n) が得られる状況では母数 1/µ の推定量は標本平均 X̄n を利用するのが一般的であるから母数 µ の推定は Yn = 1/X̄n とするのが自然である。この場合には何か問√
題はないのだろうか？標本数が大きいとき δ− 法により n[X̄n − (1/µ)] には CLT が適用できるので 1/X̄n の漸近正規性も成り立つことが期待できる。ところが X の分布が例え
15Mandelbrot, B. (1963),”The variation of some speculative prices,” Journal of Business, 26,394-419.

15

ば正規分布などであれば、推定量の評価基準としては標準的な MSE(平均二乗誤差) 基準
を用いると E[|1/X̄n |] = ∞ であるから MSE(1/X̄n ) は発散してしまうことが分かる。
この例のような状況が生じる応用例としては多変量解析、econometrics における固有値、固有ベクトルを利用する構造方程式 (structural equation) の推定量の小標本問題などが知られている16 。

注 3：第一種・第二種の合流型超幾何関数の関係について
Erdélyi [2] にしたがい関数 Φ(a, c; x) と Ψ(a, c; x) についての関係 (3.3) を補足する。関数Φ(a, c; x) は (3.1) の一つの解Φ(a, c; x) = F1 (a, c; x) = 1 +ax
a(a + 1) x+
+ ···c 1! c(c + 1) 2!

である。ここで exu を展開、ベータ関数とガンマ関数 B(n + a, c − a) = Γ(n + a)Γ(c −a)/Γ(n + c) の関係を利用すると

Γ(c)Γ(a)Γ(c − a)∫ 1
xu a−1e u(1 − u)

c−a−1

0

∫
∞
∑
Γ(c)
xn 1 n+a−1du=u(1 − u)c−a−1 duΓ(a)Γ(c − a) n=0 n! 0∞∑
xn Γ(n + a) Γ(c)=n! Γ(a) Γ(n + c)n=0∫1と表現できることに注意する。ガンマ関数を複素関数に拡張、留数の定理を Γ(x) = 0 e−t tx−1 dt+∫ ∞ −t x−1
∑∞
n
e
t
dt
の第一項
n=0 (−1) /[n!(x + n)] に適用する。ガンマ関数 Γ(−s) は −s =1
n = 0, 1, · · · に１次の極を持つので留数 (−1/2π < arg((−x) < 1/2π, 0 > γ > −Re(a)) はRes(Γ(−n)) = limx→−n,z̸=−n (x + n)Γ(x) = (−1)n /n! となるので（符号に注意して）
1 Γ(c)
2πi Γ(a)∫ γ+i∞
γ−i∞∞
Γ(−s)Γ(a + s)1
Γ(c) ∑ Γ(a + n) (−1)ns(−x) ds=
× 2πi
(−x)n
Γ(c + s)2πΓ(a) n=0 Γ(c + n) n!

Γ(c) ∑ Γ(a + n) xn
Γ(a) n=0 Γ(c + n) n!
∞=が得られ Φ(a, c; x) に一致する (Erdélyi の Page 256, (4) に対応する）。
∫∞
1
他方、Ψ(a, c; x) = Γ(a)exp(−xt)ta−1 (1 + t)c−a−1 dt より
0

Γ(a)Γ(a − c + 1)Ψ(a, c; x) =∫ ∞exp(−xt)ta−1 Γ(a − c + 1)(1 + t)c−a−1 dt

016

例えば Mariano and Sawa (1972),”The Exact Finite-Sample Distribution of the Limited-Information
Maximum Likelihood Estimator in the Case of Two Included Endogenous Variables,” Journal of theAmerican Statistical Association (JASA) が示したように固有ベクトルの推定上で生じる。次に掲示されている Kunitomo and Yuasa (2023), ”An Asymptotically Optimal Two-Sample Instrumental Variables
Estimation Using Many Instruments”の 4.1 節に関連事項が説明されている。
https://stat-expert.ism.ac.jp/wp/wp-content/uploads/2023/06/SSE-DP-2023-4.pdf

16であるが、二項展開 (1 + t)c−a−1 =

1
2πi

∑∞

n=0 (

c−a−1 n)t より
n∫ γ+i∞∞∑
1(−1)nΓ(−s)Γ(a − c + s + 1)t ds=
× 2πi
Γ(a − c + 1 − n)tn
2πi
n!
γ−i∞
n=0
s

=Γ(a − c + 1)∞∑
(a − c + 1) · (a − c − n)
n=0

n!

(−t)n

=Γ(a − c + 1)(1 + t)c−a−1となることを利用すると∫ ∞Γ(a)Γ(a − c + 1)Ψ(a, c; x)=exp(−xt)ta−1 Γ(a − c + 1)(1 + t)c−a−1 dt
0∫ γ+i∞∫ ∞1=dsΓ(−s)Γ(a − c + s + 1)
exp(−xt)t(a+s)−1 dt
2πi γ−i∞0∫ γ+i∞
1
Γ(−s)Γ(a + s)Γ(a − c + s − 1)x−s−a ds=
2πi γ−i∞∫ γ+i∞′1′′′=
Γ(−s )Γ(a + s )Γ(1 − c − s )xs ds
2πi γ−i∞′が得られる。(Erdélyi の Page 256, (5) に対応、ただし変数を a + s = −s に変更した。)
また変数を以下のように c − 1 − s = −u とすると
∫ γ+i∞1
Γ(a − c + 1 + s)Γ(−s)Γ(1 − (2 − c) − s) 1−c+s
1−cx Ψ(a − c + 1, 2 − c; x)=xds
2πi γ−i∞ Γ(a − c + 1)Γ(a − c + 1 − (2 − c) + 1)
∫ γ+i∞
Γ(−u)Γ(a + u)Γ(1 − c − u) u
1=x ds
2πi γ−i∞Γ(a)Γ(a − c + 1)∞∑
Γ(a + n)Γ(1 − c − n) (−x)n=Γ(a)Γ(a − c + 1)n!
n=0
より x1−c Ψ(a − c + 1, 2 − c; x) = Ψ(a, c; x) となる (Erdélyi の Page 257, (6) に対応)。
′′
ここで Γ(−s )Γ(1 − c − s ) の全ての極が（1 − c は整数でないとき）１次 (simple) であるので、Ψ(a, c : x)(Erdélyi の Page 256, (5) に対応) の留数を計算すると次のように Ψ(a, c; x)′′と Φ(a, c; x) の関係 Erdélyi の Page 257, (7) (3 節の (3.3)) が得られる。まず極 Γ(−s ) (s =n = 0, 1, · · · ) について留数を評価すると∞∑(−1)n nΓ(a + n)Γ(1 − c − n)x
Γ(a)Γ(a − c + 1)Ψ1 (a.c : x) =
n!
n=0
より

Ψ1 (a.c : x) =∞Γ(1 − c) ∑ Γ(a + n) Γ(1 − c − n) (−x)nΓ(a − c + 1) n=0 Γ(a)
Γ(1 − c)n!

となる。ここで符号に注意すると和の項は Φ(a, c; x) に一致するので第一項が得られる。

17

′

′

次に極 Γ(1 − c − s ) (−(1 − c − s ) = n = 0, 1, · · · ) について留数を評価すると

Ψ2 (a.c : x)=x

1−c

=x1−c

∫ γ+i∞1
Γ(a + 1 − c + u)Γ(−u − 1 + c)Γ(−u)×(−x)u du2πi γ−i∞Γ(a)Γ(a − c + 1)∞∑
Γ(a + 1 − c + n) Γ(−1 − u − n) xn
n=0

Γ(a)

Γ(a − c − +1) n!

∞Γ(c − 1) (1−c) ∑ Γ(a + 1 − c + n) Γ(c − 1 − n) xn=xΓ(a)Γ(a − c + 1)
Γ(c − 1) n!
n=0=

Γ(c − 1) (1−c)
xΦ(a − c + 1, 2 − c : x)Γ(a)と表現が得られる。ここで∫ γ+i∞1Γ(２−ｃ)
Γ(−s)Γ(a − c + 1 + u)
Φ(a − c + 1, 2 − c; x)=(−x)u ds
2πi Γ(a − c + 1 + s) γ−i∞
Γ(2 − c − 1 + u)∞
∑ Γ(a − c + 1 + n) Γ(2 − c) xn=
Γ(a − c + 1) Γ(2 − c + n) n!
n=0より第二項が得られる。第一項と第二項を合わせると公式 (3.3) (Ψ1 (a.c; x) + Ψ2 (a.c; x) に対応) が得られる。

注 4：複素関数の可視化：幾つかの例
なお参考としてガンマ関数（絶対値）、指数関数・二次関数 (正則関数の例) の図を掲載しておく17 。１変数の複素関数は 2 次元から 2 次元への写像なので実部と虚部の表現により幾つか可視化の表現がある。こうした可視化により古典的な (1 変数の) 複素関数論 (complex
analysis) の内容をより直観的に理解できる可能性がある。

17

図は西が Python によりプログラムを作成、プログラム例を付論 A に掲載した。

18

図 5: ガンマ関数 (絶対値)

図 6: 指数関数 (実部)19

図 7: 指数関数（虚部）

図 8: 二次関数（実部）

20

文献[1] Chamber, J., C. Mallows and B. Stuck, (1976). A method for simulating stable randomvariables, Journal of the Amercan Statistical Association ,Vol. 71 (354), 340-344. (Correction 82(2), 497-503).
[2] Erdélyi, A. (1953). Higher Trancendental Functions, Vol.1. New York: McGraw-Hill.
[3] Fioro, C., Hajivalliliou V. A.,& Phillps P.C.B. (2010), Bimodal t-ratios, Econometric
Journal, 13, (2), 271-289.
[4] Logan, B. F., Mallows, C. L., Rice, S. O., & Shepp, L. A. (1973). Limit distributionsof self-normalized sums. The Annals of Probability, 1(5), 788-809.
[5] Nolan, J. (2020). Univariate Stable Distributions: Models for Heavy Tailed Data,
Springer.
[6] 統計学基礎 (改訂版), 日本統計学会編, 2015, 東京図書.
[7] 現代数理統計学 , 竹村彰通, 2020, 学術図書出版.

付論 A：計算プログラム

21

# By Miwa 2024-4-30##set.seed(1)n <- 10
N <- 100000
t.cauchy <- numeric(N)for(i in 1:N) {x <- rcauchy(n)t.cauchy[i] <- sqrt(n)*mean(x)/sqrt(var(x))}#ヒストグラム.
hist(t.cauchy)#set.seed(1)n <- 30
N <- 100000
t.rt <- numeric(N)for(i in 1:N) {x <- rt(n,2)
t.rt[i] <- sqrt(n)*mean(x)/sqrt(var(x))}#ヒストグラム.
hist(t.rt)

# Nishi 2024-5-7#%%
import numpy as np
import matplotlib.pyplot as plt
from scipy import stats# %%def gen_t_stats(n, rng, df=1, n_sim=100000):
if df == np.infty:x = rng.standard_normal(size=(n, n_sim))else:
x = rng.standard_t(df=df, size=(n, n_sim))t = np.sqrt(n) * np.mean(x, axis=0) / np.sqrt(np.var(x, axis=0))
return t
rng = np.random.default_rng(123)dfs = [1, 2, 10, np.inf]
ns = [10, 30, 100, 1000]
fig, axes = plt.subplots(len(ns), len(dfs), figsize=(10, 10))for i, n in enumerate(ns):for j, df in enumerate(dfs):
axes[i, j].set_title(f"{df=}, {n=}")axes[i, j].hist(gen_t_stats(n, rng, df), bins=50)
fig.tight_layout()
fig.show()# %%rng = np.random.default_rng(123)dfs = [1.2, 1.4, 1.6, 1.8]
ns = [10, 30, 100]
fig, axes = plt.subplots(len(ns), len(dfs), figsize=(10, 10))for i, n in enumerate(ns):for j, df in enumerate(dfs):
axes[i, j].set_title(f"{df=}, {n=}")axes[i, j].hist(gen_t_stats(n, rng, df), bins=50)
fig.tight_layout()
fig.show()# %%def gen_t_stats_stable(n, rng, alpha=1, n_sim=100000):
dist = stats.levy_stable(alpha=alpha, beta=0.0)
x = dist.rvs(size=(n, n_sim), random_state=rng)t = np.sqrt(n) * np.mean(x, axis=0) / np.sqrt(np.var(x, axis=0))
return t
rng = np.random.default_rng(123)alphas = [0.1, 1, 1.5, 2]
ns = [10, 30, 100]
fig, axes = plt.subplots(len(ns), len(alphas), figsize=(10, 10))

for i, n in enumerate(ns):for j, alpha in enumerate(alphas):axes[i, j].set_title(f"{alpha=}, {n=}")axes[i, j].hist(gen_t_stats_stable(n, rng, alpha), bins=50, color="orange")
fig.tight_layout()
fig.show()# %%# Nishi 2024-5-7# %%import numpy as np
import matplotlib.pyplot as plt
from scipy import stats# %%def gen_t_stats(n, rng, df=1, n_sim=100000):
if df == np.infty:x = rng.standard_normal(size=(n, n_sim))else:
x = rng.standard_t(df=df, size=(n, n_sim))t = np.sqrt(n) * np.mean(x, axis=0) / np.sqrt(np.var(x, axis=0))
return t
rng = np.random.default_rng(123)dfs = [1, 2, 10, np.inf]
ns = [10, 30, 100, 1000]
fig, axes = plt.subplots(len(ns), len(dfs), figsize=(10, 10))for i, n in enumerate(ns):for j, df in enumerate(dfs):
axes[i, j].set_title(f"{df=}, {n=}")axes[i, j].hist(gen_t_stats(n, rng, df), bins=50)
fig.tight_layout()
fig.show()# %%rng = np.random.default_rng(123)dfs = [1.2, 1.4, 1.6, 1.8]
ns = [10, 30, 100]
fig, axes = plt.subplots(len(ns), len(dfs), figsize=(10, 10))for i, n in enumerate(ns):for j, df in enumerate(dfs):
axes[i, j].set_title(f"{df=}, {n=}")axes[i, j].hist(gen_t_stats(n, rng, df), bins=50)
fig.tight_layout()
fig.show()# %%def gen_t_stats_stable(n, rng, alpha=1, n_sim=100000):
dist = stats.levy_stable(alpha=alpha, beta=0.0)
x = dist.rvs(size=(n, n_sim), random_state=rng)t = np.sqrt(n) * np.mean(x, axis=0) / np.sqrt(np.var(x, axis=0))
return t
rng = np.random.default_rng(123)alphas = [0.1, 1, 1.5, 2]
ns = [10, 30, 100]
fig, axes = plt.subplots(len(ns), len(alphas), figsize=(10, 10))for i, n in enumerate(ns):for j, alpha in enumerate(alphas):axes[i, j].set_title(f"{alpha=}, {n=}")axes[i, j].hist(gen_t_stats_stable(n, rng, alpha), bins=50, color="orange")
fig.tight_layout()
fig.show()# %%# Nishi (2024-9-17)# %%import numpy as np
import matplotlib.pyplot as plt
from scipy import special# %%reals = np.linspace(-5, 5)
imgs = np.linspace(-5, 5, endpoint=True, num=51)# 3Dプロットを行う関数
def plot_3d_complex(func, reals, imgs, title=None, view=dict(elev=30, azim=-120)):
fig, ax = plt.subplots(
1,
1,
subplot_kw={"projection": "3d"},)xs, ys = np.meshgrid(reals, imgs)
z = xs + ys * 1.0j
surface = ax.plot_surface(xs, ys, func(z), cmap="jet", alpha=0.4)
ax.set_title(title)
ax.set_xlabel("Re")ax.set_ylabel("Im")# ax.set_zlabel(zlabel)ax.view_init(**view)fig.colorbar(surface, ax=ax, shrink=0.5)# ガンマ関数の絶対値を表示
plot_3d_complex(
lambda x: np.abs(special.gamma(x)),
reals,
imgs,
title="|Gamma|",)reals = np.linspace(0, 2)
imgs = np.linspace(-2 * np.pi, 2 * np.pi)
plot_3d_complex(
lambda x: np.real(np.exp(x)),
reals,
imgs,
title="$\Re[\exp(x)]$",
view=dict(elev=30, azim=-120),)
reals = np.linspace(-2, 2)
imgs = np.linspace(-2, 2)
plot_3d_complex(
lambda x: np.real(x**2),
reals,
imgs,
title="$\Re[x^2]$",
view=dict(elev=30, azim=-120),)

付論 B：安定分布について182 項分布の近似に端を発した中心極限定理 (CLT) は歴史的には様々な方向に拡張が試みられた。しかし、例えば確率変数の分散 E[Z 2 ] = ∞ のときは和の分布はどうなるのであろうか？例えばコーシー分布をとると E[|Z|] = ∞ で与えられる。この場合には特性関数19は φ(t) = E[exp(itZ)] = e−|t| となる（なお i2 = −1 である)。したがって、もし独立な確∑n率変数 Zi (i = 1, · · · , n) がコーシー分布にしたがっていれば和 Sn = i=1 Zi の特性関数
は φSn (t) = e−n|t| であるから確率変数 Yn = Sn /n の特性関数は φYn (t) = e−|t| となり再び√コーシ分布にしたがう。ここで注目すべき点は基準化は Sn / n ではなく Sn /n1/α (α = 1)となることである。これに対して中心極限定理 (central limit theorems, CLT) は α = 2 の場合に対応する。
コーシー分布は期待値も発散する極端な確率分布と考えられるのでもう少し現実的に期待値は定義できるが分散 E[Z 2 ] = ∞ となる場合を考察しよう。
例 A.1 : 確率分布として P (Z > z) = P (Z < −z) (z > 0) および P (|Z| > z) = z −α (z ≥1, 0 < α < 2) を取りあげる。特性関数 φ(t) を考察すると

∫ ∞
∫ −1
α
itz α −α−1
(1 − e ) z
1 − φ(t)=(1 − eitz ) |z|−α−1 dz
dz +
2
2
1
−∞
∫ ∞
1 − cos(tz)
dz
=α
z α+1
1
∫ ∞
1 − cos(u)
α
=t α
du
uα+1
t
∫
∫∞
となる。ここで u−α−1 du は可積分なので定数 C = α 0 (1 − cos(u))/uα+1 du とすると
t → 0 のとき 1 − φ(t) ∼ C|t|α と近似できる。したがって、確率変数 Zi (i = 1, · · · , n) が∑n
互いに独立にこの分布にしたがっているとき、和 Sn = i=1 Zi の特性関数は
[
]n
Sn
t
E[exp(it 1/α )] = φ( 1/α )
n
n
[
]n
t
= 1 − (1 − φ( 1/α ))n→e−C|t|α(n → ∞)となる。こうした評価より分布の裾が厚い場合には一般に確率変数和の極限分布の特性関数は定数 C2 , C3 を利用して
[
]
∫ ∞
∫ −1
dz
dz
itz
itz
(e − 1) α+1 + C3
φZ (t) = exp
C2(e − 1) α+1
z
|z|
1
−∞
となることが期待される。
18「統計エキスパート演習 2023」第 2 章付論の再掲。
https://stat-expert.ism.ac.jp/wp/wp-content/uploads/2024/01/SSE-DP-2024-1-3.pdf
19
特性関数は実数 t に対する複素関数 φ(t) = E[exp(itZ)] で与えられ、分布関数と 1-1 対応する。コー
シー分布の特性関数を求めるには複素積分を利用する必要があるが、詳しくは Durrett (1991) など確率論の教科書を参照されたい。なお最近ではコーシー分布の特性関数の導出はネット見ることができるようであるが、複素解析の基礎知識は必要のようである。

28

定義 A.1 : 確率変数 Z が無限分解可能分布 (infinitely divisible distribution) にしたがうとは、和の分布が L[Zn1 + · · · + Znn ] = L[Z] となる互いに独立・同一分布にしたがう確率変数列 Zni (i = 1, · · · , n) が存在することである。
原点での漸近挙動を考慮した表現は次のようにまとめることができる。
定理 A.1 : 無限分解可能分布にしたがう確率変数 Z の特性関数は C0 , C1 (> 0) を定数と
して
[
]
∫
C1 2
itz
φZ (t) = exp itC0 −t + (e − 1 − itz1D )ν(dz)
(A.4)
2
∫
と表現される。ここで D = {|z| ≤ 1}, (z 2 ∧ 1)ν(dz) < ∞, ν({0}) = 0 である。(ν(·) はLevy 測度と呼ばれている。)
注意 : Levy(レビー) 測度は異なる形で表現されることがある。無限分解可能分布はジャンプ確率過程と密接に関係する20 。一般に無限分解可能分布は重要な確率分布を含んでいる。例で利用した確率分布の特性関数は安定分布 (stable distribution) の特性関数の形となるが、例では中心部分を無視したことに注意する必要がある。分散が有限でない場合には分布の中心部分の評価と裾の評価を別々に行い、最終的に組み合わせる必要が生じる。
安定分布は次のよう定義され、特性関数で表現されることが知られているが詳細は省略する。
定義 A.2 : 独立で同一分布にしたがう確率変数列 Zi (i, · · · , n) が確率変数 Z の分布と同一であり、ある定数 an , bn が存在し和の分布 L[Z1 + · · · + Zn ] = L[an Z + bn ] となるとき
Z の分布は安定分布 (stable distribution) と呼ばれる。
定理 A.2 : (i) 安定分布にしたがう確率変数 Z の特性関数は C1 , C2 , C3 , 0 < α ≤ 2 を定数として次の形で与えられる。

[

∫ ∞
C1 t2
dz
φZ (t)=exp itC0 −+ C2(eitz − 1 − itz1D ) α+1
2
z
0
]
∫ 0
dz
+C3
(eitz − 1 − itz1D ) α+1|z|−∞
ただし D = {|z| ≤ 1} である。
(ii) 特性関数は C0 , d (d > 0), θ (−1 ≤ θ ≤ 1) を定数として次の形になる。
(i) 0 < α < 1, 1 < α ≤ 2 のとき
[
(
)]
t
π
α
,
φZ (t) = exp itC0 − d|t| 1 + iθ tan α
|t|
2

(A.5)

(A.6)

(ii) α = 1 のとき
[

(
)]
t 2
φZ (t) = exp itC0 − d|t| 1 + iθ
log |t|.
|t| π20

例えば佐藤健一「加法過程」(岩波書店) が詳しい。

29

(A.7)

付論 C：複素積分について 21
数理統計学や確率論では特性関数を利用すると様々な議論が一般的かつより簡明に説明できることが多い。しかしながら例えば大学文系では複素数や複素関数はあまり登場しない。
ここでは複素数と複素数の関数の有用性についての比較的分かり易い読み物として、志賀浩二「複素数 30 講」(朝倉書店) を挙げておく。複素関数を巡る議論は統計学の範囲内
でも統計的時系列分析22 (statistical time series analysis) と呼ばれる分野では有用である。
そこでは時間領域 (time domain) ではなく周波数領域 (frequency domain) での統計分析、
例えば季節性 (seasonality) を巡る時系列の周期性・循環性などの分析に役立つ。
特性関数を正確に求めるには複素数値をとる関数の積分を評価する必要がある。ここで複素数関数 f (z) の積分は実数軸と虚数軸を二次元ととらえて経路 C 上で (線積分と呼ば∫れている積分値) C f (z)dz を意味するが、通常の積分と同様に経路 C に沿って複素数 z とz の関数 f (z)(複素数) を有限和により∫
f (z)dz ∼(2B.1)Cn∑f (ξi )(zi − zi−1 )(ξi ∈ [zi−1 , zi ], i = 1, · · · , n)

i=1

して n についての極限操作により定められる。一見すると実関数の定積分とは結びつかないような複素積分の議論により、他の方法では評価が困難な積分が求められる。何回でも微分可能である性質の良い関数 (正則関数, 解析関数) であれば∫
(2B.2)
f (z)dz = 0
Cとなることを主張するコーシーの積分定理が基本的に重要である。
標準正規分布の特性関数の評価では積分経路 C : [−M, 0] → [M, 0] → [M, −it] →2[−M, −it] → [−M, 0] をとる。(図 5 を参照。) f (z) = e−z は正則関数 (regular) なので
Cauchy の積分定理を利用して閉経路 C を分割して C1 : (−M, 0) → (M, 0)(実軸上で積分値は知られている), C2 : (M, 0) → (M, −it), C3 : (M, −it) → (−M, −it)(求めたい積分の経路) C4 : (−M, −it) → (−M, 0) と云う４つの部分に積分を分割する。
このとき M → ∞ とすると (f (z) の絶対値はゼロに収束するので M を十分に大きくす
るとゼロに収束する C2 , C4 上の積分を除き）求める積分値は通常の実数軸上の積分値に帰着される。経路の方向に注意すると (C2 上の符号を変更する) コーシーの積分定理より
結局は
]
]
[
[
∫ ∞
∫ ∞
1
1
1 2
1
2
√
exp − x dx+(−1) √exp − [x − it] dx
2
2
2π −∞
2π −∞
=0

(ただし上辺での x は実数) より左辺の第一項が 1 であるから結果が得られる。
次に Cauchy 分布の特性関数の場合はより複雑である。密度関数は f (x) = 1/[π(1 + x2 )]であり、eitx = cos(tx) + i sin(tx) より sin(tx) 関数は奇関数なので積分するとゼロとなる
21「応用をめざす数理統計学」（国友直人, 朝倉) の HP 資料をほんの少し修正したものを掲示。
例えば古典的な文献として Anderson, T.W. (1971), The Statistical Analysis of Time Series, Wileyを挙げておこう。
2230

図 9: 複素平面上の経路
ので

∫ ∞
e

itx 1

1
dx = 2π 1 + x2

−∞

∫ ∞
0

1 cos txdx
π 1 + x2を評価すればよい。まず t > 0 の場合を考える。複素関数 f (z) = ez /(1 + z 2 ) は解析的ではなく点 z = i, −i において関数値が発散する。これらの２点は極と呼ばれるが、この場合にはまず t > 0 のとき、z = i における留数 (residue) Res(eitz /(1 + z 2 ), i) = e−t /(2i) と呼ばれる量 (非正則部分を部分分数分解した 1/(z − i) の係数 (1/(2i))eit(i) に対応) を利用する必要が生じる23 。すなわち t > 0 のとき 1 次の極 i についての評価により∫ ∞e−t
cos txdx=
2πi= πe−t(2B.3)
22
1+x2i
0となる24 。さらに t < 0 の場合には積分経路を点 (0, −i) を含むようにとり類似の評価を∫ ∞ tx
行うことにより 2 0 cos
dx = πet となる。t > 0 及び t < 0 の議論の両者をまとめると
1+x2
Cauchy 分布の特性関数は

ψ(t) = e−|t|(2A.4)となる。こうして導かれた特性関数の関数形より原点で微分可能でないことがわかり、期待値が存在しないことに対応する。
こうした定積分の計算の基礎となる複素関数論は美しい理論といわれている。複素数まで考察を広げるとバラバラと思われていた三角関数や指数関数、あるいは多様な関数の性
質がより深く理解でき、一見すると評価が困難な実数軸上の定積分が統一的に求めることが可能となることによるのだろう。古典的な優れた教科書として、高木貞治「解析概論」(岩波) がある。
この場合は 1 次の極 i について limz→i [(z − i)eitz /(1 + z 2 ) より求められる。一般に留数の定理 (residue
theorem) は複素関数 f が有限個の孤立特異点 ak を除き正則とするとき、ak を内部に含む閉経路 C について∫n∑f (z)dz = 2πi
Res(f, ak )
23

C

k=1

と表現される。ここで Res(f, ak ) は特異点 z = ak における留数を意味する。(例えば杉浦, Page 304 を参照。)
24
t > 0 の時は積分経路 C1 : [−R, 0] → [R, 0], C2 : [R, 0] → [0, Ri] → [−R, 0] ＞０ととると、経路内の特異点は 1 点 (0, i) である。C2 上の積分は R → ∞ のときゼロに収束、特異点の留数を計算すると留数の定理より結果が得られる。この評価を明示的に導出しているのは (数理系では常識かもしれないが) 例えば杉浦光夫「解析入門 II」東京大学出版会, IX-9 節 Page 321 などがる。

31