科研費シンポジウム

データサイエンスと周辺領域の双方向的理解への挑戦

科学研究費補助金 基盤研究(A) 20H00576「大規模複雑データの理論と方法論の革新的展開」
研究代表者 : 青嶋 誠 (筑波大学)
開催責任者 : 片山 翔太 (慶應義塾大学)
開催日時 : 2022 年 9 月 9 日(金)，9 月 10 日(土)開催場所 : 慶應義塾大学日吉キャンパス来往舎中会議室 + Zoom (ハイブリッド開催)

内容・目的 : 今やデータサイエンスはその周辺領域ー例えば因果推論，計量経済，AI・機械学習，バイオインフォマティクス，情報理論などーにおける共通言語となっており，各領域でそれぞれ独自の発展も見せています．本シンポジウムでは，統計科学の理論・応用に限らず，周辺領域の研究報告や，最先端の研究動向・問題提起などをデータサイエンスに関連する限り広く募集します．各領域における独自的な発想や，共通項について共有・議論しながら，データサイエンスをさらに深く理解し，それを周辺領域へと還元することを目指します．若手からベテランまで幅広い層の講演・参加を歓迎します．

1

プログラム9 月 9 日 (金)13:10 ~ 13:15 Opening
13:20 ~ 13:55

対面

森川 耕輔 (大阪大学), Jae Kwang Kim (アイオワ州立大学)

Semiparametric adaptive estimation under informative sampling

13:55 ~ 14:30

対面

米倉 頌人 (千葉大学), 菅澤 翔之助 (東京大学)

Adaptation of the Tuning Parameter in General Bayesian Inference with
Robust Divergence

14:30 ~ 15:05

対面高橋 慎 (法政大学), 山内 雄太 (名古屋大学), 渡部 敏明 (一橋大学),

大森 裕浩 (東京大学)
Realized stochastic volatility models with skew-t distributions

15:05 ~ 15:25 休憩
15:25 ~ 16:00

対面鈴木 洋一 (慶應義塾大学量子コンピューティングセンター)

ノイズのある量子コンピュータにおける量子振幅推定法

16:00 ~ 16:35

対面

藤森 洸 (信州大学)

The Lasso-based principal component analysis for high-dimensional
stationary time series

16:35 ~ 17:10

対面

星野 匡郎 (早稲田大学), 柳 貴英 (京都大学)

A Randomization Test for the Specification of Interference Structure

17:10 ~ 17:30 休憩
特別講演
17:30 ~ 18:30

Zoom

大津 泰介 (London School of Economics)

Likelihood inference under alternative asymptotics

2

9 月 10 日 (土)10:00 ~ 10:35

Zoom

佐々木 俊一 (岩手大学), 川崎 秀二 (岩手大学)

Erdos-Renyi ランダムグラフの大数法則に関連する漸近評価について

10:35 ~ 11:10

対面

栗栖 大輔 (横浜国立大学), 深見 陸 (東京大学), 小池 祐太 (東京大学)

深層学習による時系列データの適応的推定

11:10 ~ 11:45

対面

寺田 吉壱 (大阪大学), 山本 倫生 (大阪大学)

ベクトル量子化による大規模クラスタリングの近似法とその性質

11:45 ~ 13:10 Lunch
特別講演
13:10 ~ 14:10

対面

二宮 嘉行 (統計数理研究所)

傾向スコア解析のための情報量規準および差分の差法への展開

14:10 ~ 14:30 休憩
14:30 ~ 15:05

対面

原田 和治 (東京医科大学), 藤澤 洋徳 (統計数理研究所)

外れ値に頑健な因果効果の推定量

15:05 ~ 15:40

対面

劉 言 (早稲田大学), 谷口 正信 (早稲田大学), Hernando Ombao(King Abdullah University of Science and Technology)
Statistical Inference for Local Granger Causality

15:40 ~ 16:15

対面

山本 倫生 (大阪大学)

弱い共通サポート条件下での確率的介入に基づく因果効果の推定

16:15 ~ 16:20 Closing

3Semiparametric Adapative Estimation Under
Informative SamplingେࡕେֶେֶӃ‫ڀݚֶ޻ૅج‬Պɹ ৿઒ ߞีΞΠΦϫभཱେֶ౷‫ֶܭ‬෦
Jae Kwang Kim
1. ͸͡Ίʹ
ඪຊௐࠪͰ͸ɼඪຊΛநग़ͨ͠฼ूஂ͕ඞͣ͠΋฼ूஂશମΛ୅දͤͣɼภͬͨඪຊ (baised
sampling) ͔͠ಘΒΕͳ͍ͱ͍͏໰୊͕ੜ͡Δɽ͜ͷΑ͏ͳ໰୊͸ܽଌ஋σʔλղੳ΍ҼՌਪ࿦ʹ͓͍ͯ΋ੜ͡Δ͕ɼඪຊௐࠪʹ͓͍ͯ͸ಛʹɼσʔλ͕நग़͞ΕΔ֬཰Ͱ͋Δแ‫཰ؚ֬‬ʢinclusion
probabililtyʣͷ৘ใ΋ਪଌʹ༻͍Δ͜ͱ͕Ͱ͖ΔɽຊใࠂͰ͸แ‫཰ؚ֬‬ͷ৘ใΛ༗ޮʹ༻͍ͨɼ฼ूஂ෼෍ʹର͢ΔηϛύϥϝτϦοΫ઴ۙ༗ޮਪఆྔΛఏҊ͢Δɽ
2. ઃఆͱैདྷ๏
݁Ռม਺Λ Y , આ໌ม਺Λ X = (X1 , . . . , Xp )! ͱ͢Δɽ͋Δແ‫ݶ‬฼ूஂ͔Βͷ N ‫ݸ‬ͷແ࡞ҝඪ
ຊ FN = {(xi , yi )Ni=1 } Λߟ͑Δɽ͜͜Ͱɼ࠷ऴతʹಘΒΕΔඪຊ͸แ‫ ཰ؚ֬‬πi (i = 1, . . . , n) Ͱநग़͞Εͨ n ‫ݸ‬ͷඪຊ {(xi , yi )ni=1 } ͱ͢Δɽͨͩ͠ɼ͜ͷแ‫ ཰ؚ֬‬πi ͸આ໌ม਺΍݁Ռม਺ɼ·ͨ‫؍‬ଌ͞Εͳ͍ม਺ʹґଘ͍ͯͯ͠΋ྑ͍ (informative sampling) ͱ͠ɼN ͸ະ஌ͱ‫ط‬஌ͷ྆ํͷ৔߹Λߟ͑Δ͕ɼ౰໘ N ͸ະ஌ͱ͢Δɽ·ͨઆ໌ͷ౎߹্ɼ֬཰ม਺ Wi := πi−1 (i = 1, . . . , n)Λఆٛ͢Δɽ͜͜Ͱਪఆର৅ͱͯ࣍͠ͷ 3 ͭΛߟ͑Δ:(a) ݁Ռม਺ͷ‫ظ‬଴஋: θ = E(Y );(b) ճ‫਺܎ؼ‬: µ(x; θ) = E(Y | X = x; θ);(c) ৚݅෇͖ີ౓: f (y | x; θ).
ྫ͑͹ɼඪຊ͔Β‫ࢉܭ‬Մೳͳ (c) ʹର͢Δਪఆྔͱͯ͠ɼ࣍ͷ Hortitz-Thompson ‫ܕ‬ͷਪఆྔ
(Horvitz and Thmopson, 1952) ͕ߟ͑ΒΕΔ:N!

δi Wi Sθ (Xi , Yi ; θ) = 0.

i=1͜͜Ͱɼδi ͸༗‫ݶ‬฼ूஂ FN ͔Βநग़͞Ε͍ͯΕ͹ 1ɼͦ͏Ͱͳ͚Ε͹ 0 ΛͱΔ֬཰ม਺Ͱ͋ΓɼSθ (x, y; θ) = ∂ log f (y | x; θ)/∂θ ͸ θ ͷείΞؔ਺Ͱ͋Δɽ͜ͷਪఆํఔࣜ͸ෆภਪఆํఔࣜͱͳ͍ͬͯΔͨΊɼͦͷղ͸ θ ʹର͢ΔҰகੑΛ͕࣋ͭɼ༗ޮਪఆྔͰ͸ͳ͍ɽ

3. ηϛύϥϝτϦοΫϞσϧͱ࠷దͳਪఆํఔࣜ
ਪఆର৅ (a)–(c) ʹରͯ͠ɼ‫ڵ‬ຯͷ͋ΔύϥϝʔλҎ֎ͷಛఆΛඞཁͱ͠ͳ͍࣍ͷηϛύϥϝτ
ϦοΫϞσϧΛߟ͑Δ:f (w, x, y | δ = 1; θ, η1 , η2 , η3 ) = "="P (δ = 1 | w, x, y)f (w | x, y; η1 )f (y | x; η2 , θ)f (x; η3 )
P (δ = 1 | w, x, y)f (w | x, y; η1 )f (y | x; η2 , θ)f (x; η3 )dwdxdy
w−1 f (w | x, y; η1 )f (y | x; η2 , θ)f (x; η3 ).
w−1 f (w | x, y; η1 )f (y | x; η2 , θ)f (x; η3 )dwdxdy

ͨͩ͠ɼη1 , η2 , η3 ͸ͦΕͧΕ [w | x, y], [y | x], [x] ͷ֬཰෼෍Λ‫ن‬ఆ͢Δະ஌ͷແ‫ݩ࣍ݶ‬ύϥϝʔλͰ͋Δɽ͜͜Ͱਪఆର৅͕ (c) ͷ৔߹͸ f (y | x; η2 , θ) = f (y | x; θ) ͱͳΔ͜ͱʹ஫ҙ͢Δɽ

4

ຊ‫Ͱڀݚ‬͸ɼຊઃఆԼʹ͓͚ΔηϛύϥϝτϦοΫ઴ۙ༗ޮԼ‫ʹݶ‬౸ୡ͢ΔηϛύϥϝτϦο
Ϋ઴ۙ༗ޮਪఆྔΛಋग़ͨ͠ɽηϛύϥϝτϦοΫ༗ޮԼ‫ͱݶ‬͸ɼηϛύϥϝτϦοΫϞσϧʹ
͓͚Δਪఆྔͷ઴ۙ෼ࢄͷԼ‫͋Ͱݶ‬ΓύϥϝτϦοΫϞσϧʹ͓͚Δ Cramér-Rao ͷԼ‫ݶ‬ͷΑ͏ͳ΋ͷͰ͋ΔɽN ͕ະ஌ͷ৔߹ɼ࣍ͷਪఆํఔࣜʹ‫͍ͨͮج‬ਪఆྔ͸ηϛύϥϝτϦοΫ઴ۙ༗ޮਪఆྔͱͳΔ:

Seff =

N!

∗
δi Wi Deff
(Xi , Yi ) = 0.

(∗)

i=1

∗ ͸ਪఆର৅͝ͱʹҎԼͷ‫Ͱܗ‬ද͞ΕΔɽ
͜͜ͰɼDeff∗(a) Deff(X, Y ) = θ − Y ;1∂µ(X; θ);2
E(W ε | X) ∂θ#$
E{π̄(X, Y )Sθ (X, Y ) | X}∗(c) Deff (X, Y ) = π̄(X, Y ) Sθ (X, Y ) −.
E{π̄(X, Y ) | X}∗(b) Deff(X, Y ) =ͨͩ͠ɼε = Y − µ(X; θ), π̄(x, y) = 1/E(W | x, y) Ͱ͋Δɽ͜͜Ͱɼਪఆର৅ (b) ʹର͢Δਪఆ
ؔ਺͸ɼKim and Skinner (2013) Ͱಋग़͞Ε͍ͯΔ΋ͷͱಉҰͰ͋Δ͜ͱʹ஫ҙ͢Δɽ
·ͨɼ฼ूஂαΠζ N ͕‫ط‬஌Ͱ͋Δ৔߹ͷ࠷దͳਪఆํఔࣜ͸ɼ࣍ͷ‫Ͱܗ‬ද͞ΕΔ:

Seff =

N %
!
i=1

∗
δi Wi D̃eff(Xi , Yi ) + (1 − δi Wi )c̃∗eff

&

= 0.

(∗∗)

∗ ͱ c̃∗ ͷৄࡉʹ͍ͭͯ͸౰೔ใࠂ͢Δɽ฼ूஂαΠζ N ͕‫ط‬஌Ͱ͋Δ৔߹ɼͦͷ৘ใΛ༻͍
D̃eff
eff
ͨ (∗∗) ͷղΛ༻͍Δ͜ͱͰɼ(∗) ͷղͰ͋Δ Kim and Skinner (2013) ͷਪఆྔΛ༏ӽ͢Δਪఆྔ
͕ߏ੒ՄೳͱͳΔɽ

4. దԠతਪఆྔ
࠷దͳਪఆํఔࣜ (∗) ͸ະ஌ؔ਺Λ‫ؚ‬ΉͨΊɼ࡞‫༻ۀ‬ϞσϧΛཁ͢Δɽྫ͑͹ɼਪఆର৅͕ (c)ͷ৔߹ɼπ̄(x, y) = 1/E(W | x, y) ͱ͍͏ؔ਺ΛϞσϦϯά͠ਪఆ͢Δඞཁ͕͋Δɽͦ͜ͰɼFerrari
and Cribari-Neto (2004) ͷϕʔλճ‫ؼ‬ͷΞΠσΟΞΛར༻ͨ͠ॊೈͳϞσϦϯάํ๏ΛఏҊ͢Δɽ
ະ஌ؔ਺Λ࡞‫༻ۀ‬ϞσϧΛ༻͍ͯਪఆͨ͠΋ͷΛ༻͍ͨ࠷దͳਪఆํఔࣜ (∗) ͓Αͼ (∗∗) ͷ θʹର͢Δղ͸ɼ࡞‫༻ۀ‬Ϟσϧ͕ਖ਼͍͠৔߹ɼηϛύϥϝτϦοΫ઴ۙ༗ޮਪఆྔͱͳΓɼԾʹ‫ͬޡ‬͍ͯͨͱͯ͠΋Ұகੑͷ͋ΔਪఆྔͱͳΔɽ౰೔͸ɼ1999 ೥ʹ࣮ࢪ͞Εͨ Canadian Workplaceand Employee Survey σʔλ (Fuller, 2009) ʹఏҊख๏Λద༻͠ɼHorvitz-Thompson ‫ܕ‬ͷਪఆྔͱൺֱ͢Δ͜ͱͰɼͦͷ༗༻ੑΛࣔ͢ɽ
5. ࢀߟจ‫ݙ‬
Ferrari, S. and Chibari-Neto, F. (2004). Beta regression for modelling rates and proportions.
Journal of Applied Statistics, 31, 407–419.
Fuller, W. A. (2009). Sampling Statistics. Hoboken: John Wiley & Sons, Inc.
Horvitz, D. G. and Thompson, D. J. (1952). A generalization of sampling without replacement froma finite universe. Journal of the American Statistical Association, 47, 663–685.
Kim, J. K. and Skinner, C. J. (2013). Weighting in Survey Analysis Under Informative Sampling.
Biometrika, 100, 385–398.
Morikawa, K. and Kim, J. K. (2022). Semiparametric adaptive estimation under informative sampling.
arXiv: 2208.06039.

5

Adaptation of the Tuning Parameter in General
Bayesian Inference with Robust Divergence
Shouto Yonekura1,3 and Shonosuke Sugasawa2,3
1
2

Graduate School of Social Sciences, Chiba University

Center for Spatial Information Science, The University of Tokyo
3

Nospare Inc.

August 28, 2022

We introduce a novel methodology for robust Bayesian estimation with robust divergence (e.g., density power divergence or -divergence), indexed by tuning parameters. It
is well known that the posterior density induced by robust divergence gives highly robust
estimators against outliers if the tuning parameter is appropriately and carefully chosen.
In a Bayesian framework, one way to find the optimal tuning parameter would be using
evidence (marginal likelihood). However, we theoretically and numerically illustrate that
evidence induced by the density power divergence does not work to select the optimal tuning parameter since robust divergence is not regarded as a statistical model. To overcome
the problems, we treat the exponential of robust divergence as an unnormalisable statistical model, and we estimate the tuning parameter by minimising the Hyvarinen score.
We also provide adaptive computational methods based on sequential Monte Carlo (SMC)
samplers, enabling us to obtain the optimal tuning parameter and samples from posterior
distributions simultaneously. The empirical performance of the proposed method throughsimulations and an application to real data are also provided.
The talk is organised as follows. We first set up the framework and then show theoretically and numerically that evidence induced by density power divergence to select thetuning parameter. Instead, we propose to estimate it based on the H-score (Hyvärinen,
2005; Dawid et al., 2015) and characterise its asymptotic behaviour. As mentioned earlier,6

our method involves functions for which it is difficult to obtain an analytic representation. Therefore, we develop an adaptive and efficient Markov chain Monte Carlo (MCMC)
algorithm based on SMC samplers (Del Moral et al., 2006). Numerical applications of
the developed methodology are also provided, then conclusions and directions for future
research are stated.

References
Dawid, A. P., Musio, M., et al. (2015). Bayesian model selection based on proper scoring
rules. Bayesian analysis, 10(2):479–499. 1
Del Moral, P., Doucet, A., and Jasra, A. (2006). Sequential monte carlo samplers. Journalof the Royal Statistical Society: Series B (Statistical Methodology), 68(3):411–436. 2
Hyvärinen, A. (2005). Estimation of non-normalized statistical models by score matching.
Journal of Machine Learning Research, 6(4). 1

7

Realized Stochastic Volatility Models with Skew-t Distributions
Makoto Takahashi∗1 , Yuta Yamauchi2 , Toshiaki Watanabe3 , and Yasuhiro Omori4
1

Faculty of Business Administration, Hosei University
2

3

Graduate School of Economics, Nagoya University

Center for the Promotion of Social Data Science Education and Research, Hitotsubashi
University
4

Faculty of Economics, University of Tokyo

Abstract
Forecasting volatility and quantiles of financial returns is essential to measure the financial tail risk suchas value-at-risk and expected shortfall. There are two important aspects of volatility and quantile forecasts:
the distribution of financial returns and the estimation of the volatility. Building on the traditional stochastic
volatility model, the realized stochastic volatility model incorporates realized volatility as the precise estimatorof the volatility. Using three types of skew-t distributions, the model is extended to capture the well-known
characteristics of the return distribution, namely skewness and heavy tails. In addition to the normal andStudent’s t distributions, included as the special cases of the skew-t distributions, two of them contain the skewnormal, and hence allow more flexible modeling of the return distribution. The Bayesian estimation scheme
via a Markov chain Monte Carlo method is developed and applied to major stock indices. The empirical study
using the US and Japanese stock indices data suggests that incorporating both skewness and heavy tail to daily
returns is important for volatility and quantile forecasts.

1 SummaryFinancial volatility, defined as a standard deviation or variance of asset returns, changes stochastically over timeand hence it is important to forecast time-varying volatility for financial risk management. Time-varying volatility
used to be estimated and forecasted by two classes of time-series models using financial returns: the generalized
autoregressive conditional heteroskedasticity (GARCH) methodologies of Engle (1982) and Bollerslev (1986), andthe stochastic volatility (SV) model (Taylor, 1994). These models are consistent with volatility clustering (highvolatility persistence), and have been extended to accommodate a phenomenon called volatility asymmetry orleverage effect, i.e., the negative correlation between today’s return and tomorrow’s volatility that is observed instock markets. One such extension is the exponential GARCH (EGARCH) model proposed by Nelson (1991).
In recent years, realized volatility (RV) has replaced these models for volatility estimation. Daily RV is calculated
as the sum of squared intraday returns over a day. Andersen and Benzoni (2009) and McAleer and Medeiros∗ Corresponding author. Email: m-takahashi@hosei.ac.jp8

(2008) provide RV reviews. To predict volatility using RV, we must model its dynamics. Many studies suchas Andersen et al. (2003) have documented that daily RV may follow a long-memory process; thus, they use
autoregressive fractionally integrated moving average (ARFIMA) models (see Beran (1994) for long-memory and
ARFIMA models). The heterogeneous autoregressive (HAR) model proposed by Corsi (2009) is used more widely
for the dynamics of RV. Although it is not a long-memory model, it is known to approximate a long-memory process
well.
Although ARFIMA and HAR models have been shown to significantly improve volatility forecasts relative to
GARCH and SV models, the RV is subject to the bias caused by market microstructure noise and non-trading
hours. There are some methods for mitigating the bias in RV: the multiscale estimators by Zhang et al. (2005)
and Lan Zhang (2006), the pre-averaging approach by Jacod et al. (2009), and the realized kernel (RK) method byBarndorff-Nielsen et al. (2008, 2009). See, for example, Aït-sahalia and Mykland (2009), Ubukata and Watanabe(2014), and Liu et al. (2015) for details.
To account for the bias in RV, two classes of hybrid models have been proposed. One is the realized SV(RSV) model proposed by Takahashi et al. (2009), Dobrev and Szerszen (2010), and Koopman and Scharth (2013),
while the other is the realized GARCH (RGARCH) model and realized EGARCH (REGARCH) models proposed byHansen et al. (2012) and (Hansen and Huang, 2016), respectively. In the existing literature, the volatility forecasting
abilities of these volatility models have been compared within each class, while a comprehensive comparison acrossdifferent classes of volatility models has not been conducted extensively. An exception is Takahashi et al. (2021)
which compare volatility forecasting abilities across different classes of volatility models including the SV, EGARCH,
HAR, RSV, and REGARCH models.
To measure the financial tail risk such as value-at-risk (VaR) and expected shortfall (ES), it is important toforecast not only volatility but also quantiles of financial returns. The empirical return distribution is more peaked
and has heavier tails than the normal distribution. Such a heavy tail or leptokurtosis can be partly explained by thetime-varying volatility, but the return distribution conditional on the volatility may still be leptokurtic. Moreover,
the return distribution may also be skewed. Therefore, several types of skew Student’s t distributions have beenused in the volatility models. For example, the generalized hyperbolic (GH) skew Student’s t distribution (Aas andHaff, 2006) has been used for the SV model (Nakajima and Omori, 2012; Leão et al., 2017) and the RSV model(Takahashi et al., 2016).
Along these lines, the RSV model is extended with three types of skew- t distributions and Bayesian estimation
scheme via a Markov chain Monte Carlo method is developed. The skew-t distributions include the ones proposedby Azzalini (1985) and Fernández and Steel (1998) as well as the GH skew-t distribution. In the empirical analysis
using the Dow Jones Industrial Average (DJIA) and the Nikkei 225 (N225) data, the volatility, VaR, and ES
forecasts of the RSV models as well as the SV, EGARCH, and REGARCH models are estimated. The volatility
forecasts are evaluated by mean squared error and quasi likelihood loss functions, while the VaR and ES forecasts
are jointly evaluated by the loss function proposed by Fissler and Ziegel (2016) with the specification suggested
by Patton et al. (2019). The predictive ability test (Giacomini and White, 2006) for these loss functions and themodel confidence set (Hansen et al., 2011) are implemented to compare the forecasting performance. The results
show that the RV improves the forecasting performance significantly, the RSV models are better than the othermodels, and the skew-t distributions improve the forecasting performance. These results confirm and complement
the previous studies mentioned above.

9ϊΠζͷ͋Δྔࢠίϯϐϡʔλʹ͓͚Δྔࢠৼ෯ਪఆ๏
‫ܚ‬ጯٛक़େֶྔࢠίϯϐϡʔςΟϯάηϯλʔླ໦༸Ұྔࢠίϯϐϡʔλ͸, ༷ʑͳՊֶٕज़෼໺ʹ͓͍ͯ, ‫ݹ‬యతͳ‫ͱࢉܭ‬ൺֱ͠, ޮ཰తͳ‫ࢉܭ‬ΛՄೳʹ͢Δͱ‫ظ‬଴͞Ε͍ͯΔ. ‫ࡏݱ‬, ྔࢠΞϧΰϦζϜ͕࣮ࡍͷྔࢠίϯϐϡʔλͰ࣮ߦͰ͖Δ‫ڥ؀‬͕੔͍ͭͭ͋ΔҰํ, ‫ࡏݱ‬ͷ࣮ྔࢠίϯϐϡʔλ͸͍ΘΏΔ Noisy Intermidiate-Scale Quantum(NISQʣσόΠεͰ͋Γ, ήʔτԋࢉ਺΍ར༻ՄೳͳྔࢠϏοτͷ਺ͳͲ, ࣮༻্͍͔ͭ͘ͷ੍໿
͕͋Δ. ͦ͜Ͱ, ͜ΕΒͷ੍໿Λߟྀͨ͠ΧελϜαϒϧʔνϯ͕͍͔ͭ͘ఏҊ͞Ε͍ͯΔ.
ຊߨԋͰ͸ɺԽֶ, ۚ༥, ‫ػ‬ցֶशͳͲ, ༷ʑͳԠ༻෼໺Ͱྔࢠ‫ࢉܭ‬ͷத֩Λ୲͏ྔࢠৼ෯ਪఆ
ΞϧΰϦζϜʹয఺Λ౰ͯΔ. ಛʹ, ྔࢠৼ෯ਪఆʹΑΔϞϯςΧϧϩαϯϓϦϯάͷޮ཰Խ͸,͜ΕΒͷԠ༻ͷத֩Λͳ͢΋ͷͰ͋Γ, ͦͷॏཁੑʹؑΈ, લड़ͷํ޲ੑʹԊͬͯ, NISQ σόΠεͰ࣮ߦՄೳͳ৽͍͠ྔࢠৼ෯ਪఆΞϧΰϦζϜͷ։ൃ͕ۙ೥੝ΜʹߦΘΕ͍ͯΔ.
ྔ ࢠ ৼ ෯ ਪ ఆ ͸, ͦ ͷ ໊ ͷ ௨ Γ ຊ ࣭ త ʹ ύ ϥ ϝ ʔ λ ਪ ఆ ໰ ୊ Ͱ ͋ Γ , ྔ ࢠ ঢ় ଶ|ψ! = sin θ |good! + cos θ |bad! ʹ ‫ · ؚ‬Ε Δ, ະ ஌ ύ ϥ ϝ ʔ λ θ Λ ਪ ఆ ͢ Δ ໰ ୊ Ͱ ͋ Δ.

͜͜Ͱ |good! ͱ |bad! ͸, ௚ߦ͢Δྔࢠঢ়ଶϕΫτϧͰ͋Γ, |good! ঢ়ଶ͕‫؍‬ଌ͞ΕΔ֬཰
(0)pgood (θ) ͸, |good! ঢ়ଶͷ֬཰ৼ෯ͷೋ৐, ଈͪ sin2 θ Ͱ༩͑ΒΕΔ. ैདྷͷྔࢠৼ෯ਪఆ͸, ྔࢠৼ෯૿෯ (amplitude amplification) ͱ, Ґ૬ਪఆ (phase estimation) Λ૊Έ߹Θ࣮ͤ‫͞ݱ‬ΕΔ[1]. ྔࢠৼ෯૿෯Ͱ͸, ૿෯ԋࢉࢠ G ͕, ྔࢠঢ়ଶ |ψ! = sin θ |good! + cos θ |bad! ʹ mk ճʹ࡞༻ͨ͠ࡍ, Gmk |ψ! = sin((2mk + 1)θ) |good! + cos((2mk + 1)θ) |bad! ͷΑ͏ʹ, λʔήοτͱ(k)ͳΔঢ়ଶ |good! Λ‫؍‬ଌ͢Δ֬཰Λ pgood (θ) = sin2 ((2mk + 1)θ) ͷΑ͏ʹ૿෯͢Δ͜ͱ͕Ͱ͖Δ.

͜͜Ͱ, mk ͸, k = 0, 1, 2, · · · , M ʹ͓͚Δ૿෯ԋࢉͷ࣮ߦճ਺Ͱ͋Δ. Ұํ, ྔࢠৼ෯ਪఆͷ΋͏̍ͭͷߏ੒ཁૉͰ͋ΔҐ૬ਪఆ͸, େྔͷྔࢠήʔτૢ࡞Λඞཁͱ͢ΔͨΊ, NISQ σόΠεͰͷ࣮ߦ͸‫࣮ݱ‬తͰ͸ͳ͍. ͦ͜Ͱզʑ͸, ྔࢠৼ෯૿෯ͱ, ‫ݹ‬యత౷‫ܭ‬ख๏ (࠷໬ਪఆ) Λ૊Έ߹Θͤ, ྔࢠৼ෯ਪఆΛߦ͏ํ๏ΛఏҊͨ͠ [2]. ਤ 1 ʹఏҊख๏ͷεΩʔϜΛࣔ͢. ۩ମతʹ͸, ҟͳΔճ਺ (ਤ 1 Ͱ͸, m0 , m1 , · · ·, mM ճ) ͷྔࢠৼ෯૿෯ΛߦͳͬͯಘͨΕͨ݁ՌΛ૊Έ߹Θͤ࠷໬ਪఆΛߦ͍, ύϥϝʔλਪఆΛߦ͏. ͜ͷख๏Λ༻͍ͯϞϯςΧϧϩੵ෼Λߦͬͨͱ͜Ζ, ‫ݹ‬యϞϯςΧϧϩੵ෼ͱൺֱ͠, ߴ଎ʹ࣮ߦՄೳͰ͋Δ͜ͱ (ೋ৐Ճ଎͕͋Δ͜ͱ) Λ֬ೝ͠, େ෯ʹྔࢠήʔτૢ࡞Λ࡟‫͖Ͱݮ‬Δ͜ͱΛࣔͨ͠.
ߋʹ, զʑ͸͜ͷํ๏Λ֦ு͠, ‫࣮ݱ‬తͳϊΠζޮՌΛߟྀʹೖΕ, ௒఻ಋ IBM ྔࢠσόΠεΛ༻͍࣮ͨ‫ݧ‬ͷղੳΛߦͬͨ. ঘ, ৄࡉʹ͍ͭͯ͸ߨԋͰड़΂Δ. ࣮‫݁ݧ‬Ռͷղੳʹ͋(k)ͨΓ, զʑ͸࣍ͷΑ͏ͳϞσϧΛಋೖͨ͠. |good! ঢ়ଶΛ‫؍‬ଌ͢Δ֬཰͕, pgood (θ, βk ) ‫ͼٴ‬,(k)qgood (θ, βk ) Ͱ͋Δ 2 छྨͷྔࢠճ࿏Λ༻ҙ͢Δ. ୠ͠, βk ͸ϊΠζͷޮՌΛද͢ύϥϝʔλͰ͋
(k)

(k)

Γ, pgood (θ, βk ) = 1/2 − 1/2βk cos(2(2mk + 1)θ), qgood (θ, βk ) = 1/2 − 1/2βk cos(2(2mk − 3)θ)ͱͨ͠. βk ͸, 0 ≤ β ≤ 1 ͷ࣮਺Ͱ͋Γ, βk = 1 ͕ ϊΠζͷͳ͍৔߹ʹରԠ͢Δ. զʑ

10

ਤ1

࠷໬ਪఆΛ༻͍Δྔࢠৼ෯ਪఆεΩʔϜ. ҟͳΔճ਺ (m0 , m1 , · · · , mM ) ͷྔࢠৼ෯૿

෯ (ࠨ). ͦΕͧΕͷ࣮‫݁ݧ‬Ռ͔ΒಘΒΕΔ໬౓ؔ਺ (தԝ). ҟͳΔճ਺ͷྔࢠৼ෯૿෯Λ૊Έ߹Θͤͨ໬౓ؔ਺ (ӈ).

(k)

(k)

͸, ͜ͷ pgood (θ, βk ) ‫ͼٴ‬, qgood (θ, βk ), ୠ͠ k = 0, 1, 2, · · · , M Λ༻͍ͯ໬౓ؔ਺Λߏ੒͠,࠷໬ਪఆΛߦͬͨ. ͜͜Ͱ໰୊ͱͳΔͷ͸, ֬཰෼෍Ϟσϧʹଟ͘ͷϊΠζύϥϝʔλ βk ,k = 0, 1, 2, · · · , M (‫֎ہ‬ύϥϝʔλ) ͕‫·ؚ‬ΕΔͨΊ, ࠷దͳύϥϝʔλͷਪఆ͕ޮ཰తͰͳ

͍ͱ͍͏ຊ࣭తͳ໰୊͕͋Δ. ͜ͷ໰୊ʹରͯ͠զʑ͸, ύϥϝʔλ௚ަԽ๏ (parameterorthogonalization method)[3, 4] Λద༻͠, ଞͷϊΠζύϥϝʔλ βk , k = 0, 1, 2, · · · , M , Λআ‫ͯ͠ڈ‬, ໨తͷৼ෯ύϥϝʔλ θ ͷΈͷਪఆΛߦͬͨ [5]. ߨԋͰ͸, ຊख๏ͷৄࡉ΍, ࣮ࡍͷ௒఻ಋྔࢠσόΠεΛ༻͍࣮ͨ‫ݧ‬ͷղੳ݁ՌͳͲʹ͍ͭͯ঺հ͢Δ༧ఆͰ͋Δ.

ࢀߟจ‫ݙ‬[1] G. Brassard, P. Høyer, M. Mosca, A. Tapp, Quantum amplitude amplification and estimation. Contemp. Math. Ser. Millenn. 305, 53-74 (2002)[2] Y. Suzuki, S. Uno, R. Raymond, T. Tanaka, T. Onodera, and N. Yamamoto, Amplitude
estimation without phase estimation, Quant. Info. Proc. 19, 75 (2020).
[3] D. R. Cox and N. Reid, Parameter orthogonality and approximate conditional inference,
J. R. Stat. Soc., Ser. B 49(1), 1 (1987).
[4] J. Suzuki, Nuisance parameter problem in quantum estimation theory: Tradeoff relationand qubit examples, J. Phys. A: Math. Theor. 53, 264001 (2020).
[5] T. Tanaka, S. Uno, T. Onodera, N. Yamamoto, and Y. Suzuki, Noisy quantum amplitude
estimation without noise estimation, Phys. Rev. A 105, 012411 (2022).

11

The Lasso-based principal component analysis forhigh-dimensional stationary time series
藤森

信州大学

1

洸

経法学部

序論

高次元の設定において，従来の主成分分析は精度が悪いことが知られており，特に，スパー
スな主成分ベクトルの推定には，より良い性質を持った推定量が，主に独立標本の設定において提案されてきた．本講演では Gaussian よりも裾の重いケースを含む高次元定常時系列に対して Lasso 型のスパース主成分分析手法を適用し，得られる推定量の漸近的な挙動について議論する．さらに，従来の主成分分析手法や独立標本に対する Lasso 型主成分分析との比較を行っていく．

2

主結果{Xt }t∈Z を確率空間 (Ω, F, P ) 上の Rp -値，平均 0 の定常過程とする．観測系列 X1 , . . . , Xn , n ∈N に対して，次の p × p 標本共分散行列及び，共分散行列を考える．
n

Σ̂n =

1!
Xt Xt" ,nΣ0 = E[Xt Xt" ].

t=1

以下，Σ0 の最大固有値 φ2max に対応する，規格化された第一主成分ベクトル q 0 の推定問題
を考える．推定の対象となるパラメータは，β 0 = φmax q 0 として特徴づけることができる
が，これは，次の最適化問題の解として与えられる：1
β 0 = arg min #Σ0 − ββ " #2F ⇐⇒ Σ0 β 0 = #β 0 #22 β 0 ,β 4
ただし，# · #F は行列のフロベニウスノルムである．さて，β 0 の推定量として，次の Lasso型推定量を定義する (van de Geer (2016))：
"
#
1
1
" 2
β̂n := arg min#Σ̂n − ββ #F + λ1 #β#1 , B := {β : #β − β 0 #2 ≤ η}.
β∈B 4

12

ただし，λ1 は正則化パラメータであり，η は適当な定数である．本講演では，適当な定数
q > 0 に対して，p = qn が成立することを仮定し，推定量 β̂n1 の挙動について議論する．
時系列における従属構造の影響は，混合係数を用いて評価することになる．例えば，ガウス型過程 {X}t∈Z に対しては次の ρ-混合係数を用いた評価を行う．

ρ(l) := sup{| Cov(f (Xt ), g(Xt+l ))| :E[f ] = E[g] = 0, E[f 2 ] = E[g 2 ] = 1},l ∈ Z.

ρ-混合係数を用いることで，ガウス型確率過程の標本共分散行列の推定誤差評価に際して，Hanson-Wright の不等式と呼ばれる集中不等式を利用することができる．いま，β 0 の非零成分の個数を s0 とする．このとき，共分散行列の固有値に対するスパイク性に関する条件や混合係数に関する条件をはじめとする適当な条件の下で，推定量の収束レートは次のように導出される．

Theorem
2.1. 確率過程 {Xt } はガウス型であるとし，その ρ-混合係数を ρ(·) とする．，
$∞
ρ(l)
<
∞, φ2max = O(1), n → ∞ を仮定し，正則化パラメータ λ1 は次を満たすと
l=0
する：
%
log p
λ1 * s 0.
n
このとき，適当な正則条件の下で，次が成立する．
&
'
%
log
p
3/2
#β̂n1 − β 0 #1 = Op s0,n#β̂n1 − β 0 #22 = Op

(

log p
s30
n

)

,

n→∞
n → ∞.

上述の定理において，ρ-混合係数の総和に関する条件を仮定しているが，これは推定量の非漸近的な評価を行った際に，時系列の従属構造の影響が現れることに由来している．すなわち，独立標本のケースとは異なった挙動を観察することができる．
本講演ではガウス型確率過程だけではなく，より裾の重いケースを含むクラスである
sub-Weibull 型確率過程をより強い混合係数に対する条件の下で扱い，Lasso 型推定量の漸近的挙動を比較する．なお，講演内容は，Fujimori et al. (2016) に基づく．

References
Fujimori, K., Goto, Y., Liu, Y., Taniguchi, M. (2021). Sparse principal component analysis
for high-dimensional stationary time series. Submitted.
van de Geer, S. A. (2016). Estimation and testing under sparsity. Springer.

13

A Randomization Test for the Specification of Interference Structure
Tadao Hoshino∗ and Takahide Yanagi†
Recent studies in the causal inference literature have highlighted the importance of accounting
for treatment spillovers from other units via empirical applications in many fields, including politicalscience, epidemiology, education, economics, etc. In a most general case, the treatment of all units
in the population may affect one’s outcome. The principle of causal inference is to compare similar
individuals with different treatment status. However, when treatment spillovers are present, since eachindividual has his/her unique social network with others, no such similar counterparts might exist.
Thus, it is generally not possible to identify meaningful causal parameters without some simplifying
assumptions on the interference structure. Indeed, Imbens and Rubin (2015) states ”causal inferenceis generally impossible without such assumptions, and thus it is critical about their content and theirjustifications.”
To address this problem, a common approach in the literature is to assume the existence ofan exposure mapping that summarizes the impacts from others’ treatments into lower dimensionalsufficient statistics. For example, Hong and Raudenbush (2006) study the impact of school retention onlater academic performance, assuming that other students’ retention may affect one’s own performance
only through whether their school has a higher or lower retention rate. As another example, Leung(2020) considers a treatment spillover model in which one’s potential outcome is a function not onlyof one’s treatment but also of the number of treated neighbors and that of all neighbors.
Clearly, if the chosen exposure mapping is not appropriate, the resulting causal inference may
be misleading (e.g., failure of detecting treatment spillovers).1 As stated in the quote from Imbens
and Rubin (2015), in order for the obtained spillover effects to be valid, it is required to justify the
specification of the exposure mapping. A direct approach to accomplish this is to statistically test the
specification, which is exactly the aim of this study. That is, we propose a randomization specification
test for the form of general exposure mappings.
We are not the first to develop this type of randomization test; see, e.g., Athey et al. (2018),
Basse et al. (2019), and Puelz et al. (2022). The major difference between our approach and theirs
lies in the construction of the test statistic. The aforementioned studies are almost agnostic about
the construction of test statistics with sufficient power, or they only provide examples that work well∗
School of Political Science and Economics, Waseda University, 1-6-1 Nishi-waseda, Shinjuku-ku, Tokyo 169-8050,
Japan. Email: thoshino@waseda.jp.
†
Graduate School of Economics, Kyoto University, Yoshida Honmachi, Sakyo, Kyoto, 606-8501, Japan. Email:yanagi@econ.kyoto-u.ac.jp.
1
Some recent studies uncover under what conditions one can estimate meaningful causal parameters even when the
exposure mapping is misspecified or is not explicitly specified (Hoshino and Yanagi, 2021; Sävje et al., 2021; Leung,2022). A common finding in these studies is that, only if the network dependence is sufficiently weak, we can identify
some composite causal parameters.

14

for certain situations. In contrast, by introducing the notion of coarseness of exposure mappings, we
propose a method that automatically derives test statistics equipped with reasonable power in mostsituations.
The basic idea is as follows. Suppose that we have two exposure mappings E 0 and E 1 . We saythat E 0 is coarser than E 1 if there exists a surjective mapping c such that c(E 1 ) = E 0 . For example,
E 0 = 1t# treated neighbors ° 0u,

E 1 = # treated neighbors

Now, suppose that we would like to test the following null hypothesis:
H0 : E 0 is a correct exposure mapping.
Then, if H0 is true, E 1 is also a correct exposure mapping. This enables us to perform a randomization
test by shuffling the treatment assignment while keeping the value of E 0 fixed but allowing the valueof E 1 to alternate. Under H0 , the distributions of potential outcomes for different E 1 -values should
be identical; otherwise, H0 should be rejected. The result of numerical simulations shows that ourapproach works satisfactorily compared to the existing method.

References
Athey, S., Eckles, D., and Imbens, G.W., 2018. Exact p-values for network interference, Journal of the American
Statistical Association, 113 (521), 230–240.
Basse, G.W., Feller, A., and Toulis, P., 2019. Randomization tests of causal effects under interference, Biometrika,
106 (2), 487–494.
Hong, G. and Raudenbush, S.W., 2006. Evaluating kindergarten retention policy: A case study of causal inference formultilevel observational data, Journal of the American Statistical Association, 101 (475), 901–910.
Hoshino, T. and Yanagi, T., 2021. Causal inference with noncompliance and unknown interference, arXiv:2108.07455.
Imbens, G.W. and Rubin, D.B., 2015. Causal Inference in Statistics, Social, and Biomedical Sciences, Cambridge
University Press.
Leung, M., 2020. Treatment and spillover effects under network interference, The Review of Economics and Statistics,
102 (2), 368–380.
Leung, M., 2022. Causal inference under approximate neighborhood interference, Econometrica, 90 (1), 267–293.
Puelz, D., Basse, G., Feller, A., and Toulis, P., 2022. A graph-theoretic approach to randomization tests of causal effects
under general interference, Journal of the Royal Statistical Society: Series B, 84 (1), 174–204.
Sävje, F., Aronow, P.M., and Hudgens, M.G., 2021. Average treatment effects in the presence of unknown interference,
The Annals of Statistics, 49 (2), 673–701.

15

LIKELIHOOD INFERENCE UNDER ALTERNATIVEASYMPTOTICSTAISUKE OTSU

In this talk, I discuss a unified inference framework for various alternative or
non-standard asymptotic inference problems, including sparse network, degeneracy
under multiway data, small bandwidth for semiparametric inference, many weak
instruments, and many covariates asymptotics. The talk is based on three papers(1) Empirical likelihood for network data (with Yukitoshi Matsushita)(2) Multiway empirical likelihood (with Harold Chiang and Yukitoshi Matsushita)(3) Jackknife empirical likelihood: small bandwidths, sparse network and highdimensional asymptotics (with Yukitoshi Matsushita)
The first paper is summarized as follows. Analysis on network data is becoming
increasingly important in various fields of data science, and the literature on statistical modelling and estimation algorithms for networks is rapidly growing. However,general statistical inference methods for networks are still less developed. This article develops concept of nonparametric likelihood for network data based on the
network moments, and proposes general inference methods by adapting the theoryof jackknife empirical likelihood. Our methodology can be used not only to conduct
inference on population network moments and parameters in network formation
models, but also to implement goodness-of-fit testing, such as testing block size forstochastic block models. Theoretically we show that the jackknife empirical likelihood statistic loses its asymptotic pivotalness under the sparse network asymptotics
and develop a modified statistic which converges to a chi-squared distribution underboth the sparse and dense network asymptotics.
The second paper is summarized as follows. This paper develops a general methodology to conduct statistical inference for observations indexed by multiple sets of entities. We propose a novel multiway empirical likelihood statistic that converges to achi-square distribution under the non-degenerate case, where corresponding Hoeﬀding type decomposition is dominated by linear terms. Our methodology is related
to the notion of jackknife empirical likelihood but the leave-out pseudo values areconstructed by leaving out columns or rows. We further develop a modified version

16

of our multiway empirical likelihood statistic, which converges to a chi-square distribution regardless of the degeneracy, and discover its desirable higher-order property
compared to the t-ratio by the conventional Eicker-White type variance estimator.
The proposed methodology is illustrated by several important statistical problems,such as bipartite network, two-stage sampling, generalized estimating equations, andthree-way observations.
The third paper is summarized as follows. This paper sheds light on inference
problems for statistical models under alternative or nonstandard asymptotic frameworks from the perspective of jackknife empirical likelihood (JEL). Examples include
small bandwidth asymptotics for semiparametric inference, many covariates asymptotics for regression models, many-weak instruments asymptotics for instrumental
variable regression, and sparse network asymptotics. We first establish Wilks’ theorem for the JEL statistic on a general semiparametric inference problem under theconventional asymptotics. We then show that the JEL statistics lose asymptotic
pivotalness under the above nonstandard asymptotic frameworks, and argue thatthese phenomena are understood as emergence of Efron and Stein’s (1981) bias of
the jackknife variance estimator in the first order. Finally we propose a modificationof JEL to recover asymptotic pivotalness under both the conventional and nonstandard asymptotics. Our modification works for all above examples and provides aunified framework to investigate nonstandard asymptotic problems.

17Erdős-Rényi ランダムグラフの大数法則に関連する漸近評価について

岩手大学総合科学研究科理工学専攻
岩手大学理工学部

佐々木俊一
川崎秀二

はじめに：Erdős-Rényi ランダムグラフ

1

ランダムグラフ（以下，RG）は n 個のノード（頂点）とそれらの間を繋ぐエッジ（辺）から成り，ノード間に確率 p ∈ (0, 1) でエッジが形成される．1 つのノードがエッジにより連結するノード数は非負整数値の確率変数で，その期待値を λ > 0 とすると，λ = np の関係が仮定される．RGは λ = 1 を境界として，次のような相転移を呈する事が知られている：
 λ < 1： 劣臨界 · · · · 確率１で有限ノード数 (n < ∞) の系となる（「絶滅」確率 η = 1）


λ = 1： 臨界 · · · · 条件によりη = 1 あるいはη = 0 の系となる

λ > 1： 優臨界 · · · · 確率 ζ = 1 − η で無限ノード数 (n → ∞) の系となる（「生存」確率 ζ > 0）

RG のうち特に，1 つのノードからエッジにより連結するノードの数が独立同分布で，それぞれPoisson 分布 Pois(λ) に従うものを Erdős-Rényi RG（ERG）と呼ぶ．本研究では，ERG で
n → ∞ の時の漸近的性質について調べる．

また一般に，RG においてエッジによって連結したノードの集合 (クラスター) が複数個生成される. 全ノード数 n の ER RG において生成された複数のクラスターのうち，ノード数が最も多い
クラスターを最大クラスターと呼び，Cmax と表す．|Cmax | は最大クラスターのノード数である．
優臨界での最大クラスターサイズについては，次の形の大数の法則が成り立つ；!
定理 1 ([1])

"

#

$

λ > 1 とする．任意の ν ∈ ( 12 , 1) に対して，δ = δ(ν, λ) が存在して，次を満たす：
%&
'
&
(
)
(1)
P & |Cmax | − ζλ n & ! nν = O n−δ ,as n → ∞.

図 2. (1) で ν を固定した時の P 対 n

図 3 (1) で ν = 0.51 を固定した時の log P 対 λ18

n が大きくなると (1) 式左辺の確率 P が小さくなるのは当然だが，図３から分かるように，λ に関しても，P は |λ − 1| の減少関数となっている事が見受けられる．この事を説明できるよう，Pを λ の関数として評価する事を目指す．

さて，上記の定理１に対応する中心極限定理 (CLT) も，ある σ 2 = σλ2 に対し，次のように知られている [1]：()|Cmax | − ζn√=⇒ N 0, σ 2 ,
n

(2)

as

n → ∞.

つまり，個々のミクロ因子である各ノードのリンク数の分布も分散が有限であり，またマクロにもCLT が成り立つような分散が有限の系である．では NW の構造やサイズは，単純に独立確率変数
列の和として捉えられるかというと一般にはそうではなく，各ノードからのリンクが現時点の最大クラスタに接続するか否かのランダム性がある．接続すれば，その新規ノードが Cmax へ追加され，|Cmax | が増大する．そのランダム性の結果，|Cmax | は優臨界では n の線形関数 ζn に漸近す
|Cmax |
る．劣臨界では，−→ Iλ−1（確率収束）となる事が知られている．ただし Iλ = λ − 1 − log λ
log n
である．

|Cmax | 自体の分布は，必ずしも正規分布に従わない．幾つかの λ に対する |Cmax | のヒストグラ
ムの変化を図 4 に示す．
1.8

10 -3

10 -3

0.16

4.5

|Cmax||Cmax||Cmax|

1.6

4

0.14

1.4

3.5

0.12

1.2

3

1

2.5

0.8

2

0.6

1.5

0.4

1

0.2

0.5

0

0

0.1
0.08
0.06

500

1000

1500

2000

|Cmax|

図 4(a): λ = 1

図４：

2500

0.04
0.02

1.56

1.58

1.6

0

1.62

|Cmax|

図 4(b): λ = 2

10 4

1.99851.999

1.9995|Cmax|2
10 4

図 4(c): λ = 8

幾つかの λ に対する P のヒストグラム．λ が小さい方から大きい方へ向かって，
右歪み → 左歪みと変化している．ビンの幅が異なるのは，横軸のレンジが変化している事による．

図４から，|Cmax | の分布は，λ が小さい方から右歪み → 左右対称 → 左歪みと変化している．
左右対称のケースでは，正規分布が良くフィットしている．また，右歪みあるいは左歪みのうち，
パラメータの値によっては対数正規分布がフィットする時もあれば，いずれもフィットしない時もある．これらの分布の変化について，劣臨界／臨界／優臨界を包摂した分布の確率的メカニズム
や構造について調べることを目標としている．

参考文献[1] R. van der Hofstad, Random Graphs and Complex Networks, Cambridge Univ. Press, 2017.

19

ADAPTIVE DEEP LEARNING FOR NONPARAMETRIC TIME SERIESREGRESSIONDAISUKE KURISU, RIKU FUKAMI, AND YUTA KOIKE
Abstract. In this paper, we develop a general theory for adaptive nonparametric estimation ofmean functions of nonstationary and nonlinear time series using deep neural networks (DNNs). We
first consider two types of DNN estimators, non-penalized and sparse-penalized DNN estimators,
and establish their generalization error bounds for general nonstationary time series. We then derive
minimax lower bounds for estimating mean functions belonging to a wide class of nonlinear autoregressive (AR) models that include nonlinear generalized additive AR, single index, and threshold
AR models. Building upon the results, we show that the sparse-penalized DNN estimator is adaptive and attains the minimax optimal rates up to a poly-logarithmic factor for many nonlinear
AR models. Through numerical simulations, we demonstrate the usefulness of the DNN methods
for estimating nonlinear AR models with intrinsic low-dimensional structures and discontinuous orrough mean functions, which is consistent with our theory.

1. Introduction
Let (⌦, G, {Gt }t 0 , P) be a filtered probability space. Consider the following nonparametric timeseries regression model:
Yt = m(Xt ) + ⌘(Xt )vt , t = 1, . . . , T,

(1.1)

where T
3, (Yt , Xt ) 2 R ⇥ Rd , and {Xt , vt }Tt=1 is a sequence of random vectors adapted to
the filtration {Gt }Tt=1 . We assume C⌘ := supx2[0,1]d |⌘(x)| < 1. In this paper we investigate
nonparametric estimation of the mean function m on the compact set [0, 1]d , that is, f0 := m1[0,1]d .
2. Deep neural networks
To estimate the mean function m of the model (1.1), we fit a deep neural network (DNN) witha nonlinear activation function : R ! R. The network architecture (L, p) consists of a positive
integer L called the number of hidden layers or depth and a width vector p = (p0 , . . . , pL+1 ) 2 NL+2 .
A DNN with network architecture (L, p) is then any function of the formf : Rp0 ! RpL+1 , x 7! f (x) = AL+1

L

AL

L 1

···

1

A1 (x),(2.1)

where A` : Rp` 1 ! Rp` is an affine linear map defined by A` (x) := W` x + b` for given p` 1 ⇥ p`weight matrix W` and a shift vector b` 2 Rp` , and ` : Rp` ! Rp` is an element-wise nonlinear
activation map defined as ` (z) := ( (z1 ), . . . , (zp` ))0 . We assume that the activation functionis C-Lipschitz for some C > 0, that is, there exists C > 0 such that | (x1 )(x2 )|  C|x1 x2 |for any x1 , x2 2 R. Examples of C-Lipschitz activation functions include the rectified linear unit
D. Kurisu is partially supported by JSPS KAKENHI Grant Number 20K13468. Y. Koike is partially supported
by JST CREST Grant Number JPMJCR2115 and JSPS KAKENHI Grant Number 19K13668.

20

(ReLU) activation function x 7! max{x, 0} and the sigmoid activation function x 7! 1/(1 + e x ).
For a neural network of the form (2.1), we define✓(f ) := (vec(W1 )0 , b01 , . . . , vec(WL+1 )0 , b0L+1 )0
where vec(W ) transforms the matrix W into the corresponding vector by concatenating the column
vectors.
3. Main results
We establish• generalization error bounds for non-penalized DNN and sparse-penalized DNN estimator,• minimax lower bounds for estimating mean functions f0 of a nonlinear AR model thatbelong to (i) composition structured functions and (ii) `0 -bounded affine class,• minimax optimality of sparse-penalized DNN estimator over those two classes of meanfunctions.
In addition to the theoretical results, we also conduct simulation studies to investigate the finite
sample performance of the DNN estimators and compare their performance with other estimators(kernel ridge regression, k-nearest neighbors, and random forest). We find that the DNN methods
work well for the models with (i) intrinsic low-dimensional structures and (ii) discontinuous orrough mean functions. These results are consistent with our main results.
References(D. Kurisu) Graduate School of International Social Sciences, Yokohama National University, 79-4
Tokiwadai, Hodogaya-ku, Yokohama 240-8501, Japan.
E-mail address: kurisu-daisuke-jr@ynu.ac.jp(R. Fukami) Graduate School of Mathematical Science, The University of Tokyo, 3-8-1 Komaba,
Meguro-ku, Tokyo 153-8914, Japan.
E-mail address: rick.h.azuma@gmail.com(Y. Koike) Graduate School of Mathematical Science, The University of Tokyo, 3-8-1 Komaba,
Meguro-ku, Tokyo 153-8914, Japan.
E-mail address: kyuta@ms.u-tokyo.ac.jp21ϕΫτϧྔࢠԽʹΑΔେ‫ن‬໛ΫϥελϦϯάͷۙࣅ๏ͱͦͷੑ࣭ࣉా٢ұ 1,3 , ࢁຊ ྙੜ 2,3
1 େࡕେֶେֶӃ‫ڀݚֶ޻ૅج‬Պ, 2 େࡕେֶେֶӃਓؒՊֶ‫ڀݚ‬Պ
3 ཧԽֶ‫৽ֵॴڀݚ‬஌ೳ౷߹‫ڀݚ‬ηϯλʔ

1. ͸͡Ίʹ
ۙ೥ɼσʔλͷେ‫ن‬໛ԽͱෳࡶԽ͕ਐΈ, σʔλ͔ΒԾઆ΍༗ӹͳ৘ใΛ֫ಘ͢Δ͜ͱ͕՝୊ͱͳ͍ͬͯΔ. ͦͯ͠, ୳ࡧతͳσʔλղੳͷதͰ΋‫ֶ͠ͳࢣڭ‬शͷॏཁੑ͕࠶ೝࣝ͞Ε͍ͯΔɽΫϥελϦϯά๏͸, σʔλͷഎ‫ޙ‬ͷΫϥελߏ଄Λ໌Β͔ʹ͢ΔͨΊͷ‫ֶ͠ͳࢣڭ‬शͷํ๏Ͱ͋Γ, ༷ʑͳ෼໺Ͱ޿͘Ԡ༻͞Ε͍ͯΔ. ࠷΋୅දతͳΫϥελϦϯά๏ͱͯ͠, k-means ๏͕‫͛ڍ‬ΒΕΔ. k-means ๏͸ͦͷ؆ศੑͱ‫ࢉܭ‬ίετͷ௿͔͞Βଟ༻͞ΕΔ͕, ͦͷ୯७͞‫ʹނ‬σʔλͷഎ‫͋ʹޙ‬ΔෳࡶͳΫϥελߏ଄Λे෼ʹଊ͑ΒΕͳ͍Մೳੑ͕͋Δ. ͦͷͨΊ, spectral clustering(von Luxburg, 2007) ͳͲ, ΑΓॊೈʹΫϥελߏ଄Λଊ͑ΒΕΔํ๏ͷར༻͕๬·͍͠. ͔͠͠,͜ΕΒͷํ๏͸, Ұൠʹ‫ࢉܭ‬ίετ͕ߴ͘, ‫ࢉܭ‬ίετͷ࡟‫͕ݮ‬େ͖ͳ՝୊ͱͳ͍ͬͯΔ.
Χʔωϧ๏ʹ͓͚ΔΫϥελϦϯά๏ʹରͯ͠͸, Χʔωϧ๏ʹಛԽͨۙ͠ࣅ๏Ͱ͋Δ Nyströmۙࣅ΍ Random Fourier Feature ͳͲͷ‫ࢉܭ‬ίετ؇࿨๏͕ద༻ग़དྷΔ. ҰํͰ, ΫϥελϦϯά๏ʹ͓͍ͯ, ͜ΕΒͷํ๏ʹΑΔۙࣅ͸҆ఆతʹ‫ػ‬ೳ͠ͳ͍ͱ͍͏໰୊఺͕͋Δ. Yan et al. (2009)Ͱ͸, spectral clustering ͷ k-means ๏ʹ‫ࣅۙͮ͘ج‬๏ (KASP) ͕ఏҊ͞Ε͍ͯΔ. KASP ͸,
spectral clustering ʹ‫ݶ‬Βͣ೚ҙͷΫϥελϦϯά๏ʹରͯ͠ద༻͢Δ͜ͱ͕Ͱ͖Δ. ͜ͷํ๏ͷେ͖ͳར఺͸, ҆ఆੑ, ؆ศੑ, ‫ࢉܭ‬ίετͷ௿͞, ൚༻ੑͰ͋Δ.
ຊൃදͰ͸, KASP ͷ໰୊఺Λ໌Β͔ʹ͠, ͦͷ໰୊఺Λ؇࿨ͨ͠ҰൠͷϕΫτϧྔࢠԽ (vector
quantization) Λ༻͍ͨۙࣅํ๏ΛఏҊ͢Δ.
2. KASP ͷ໰୊఺ͱఏҊख๏
KASP Ͱ͸, Ϋϥελ਺Λଟ͘ઃఆͨ͠ k-means ๏Λେ‫ن‬໛σʔλʹద༻͠, ಘΒΕͨΫϥελத৺ΛσʔλΛ୅ද͢Δ఺ͱ͢Δ. ͦͯ͠, ಘΒΕͨ୅ද఺ͷΈʹ spectral clustering ͳͲෳࡶͳ
ΫϥελϦϯά๏Λద༻͠, ୅ද఺ʹର͢ΔϥϕϧΛͦͷ୅ද఺ʹ͍ۙσʔλ఺ͷϥϕϧͱ͢Δ.
͔͠͠, k-means ๏ʹΑͬͯੜ੒ͨ͠୅ද఺͔Βߏ੒͞ΕΔ‫ݧܦ‬෼෍͸, ฼ूஂ෼෍Λ୅ද͢Δ఺ͱ͸ͳΒͳ͍. ࣮ࡍʹ, ୅ද఺ͷ‫ݧܦ‬෼෍͸, ฼ूஂ෼෍ΑΓ΋੄͕ॏ͍෼෍ʹऩଋ͢Δ͜ͱ͕ࣔ
ͤΔ. ͦͷͨΊ, KASP Λ༻͍ͯۙࣅΛߦ͏ͱ, άϥελϦϯά݁ՌʹζϨ͕ੜͯ͡͠·͏.
͜ͷ໰୊఺ʹର͢Δ࠷΋γϯϓϧͳղܾ๏͸, ֤୅ද఺ʹద੾ͳॏΈΛ༩͑Δ͜ͱͰ͋Δ. Kmeans ๏ʹରԠ͢ΔϕΫτϧྔࢠԽ͸, ฼ूஂ෼෍ͱͷ L2 -Wasserstein ‫཭ڑ‬Λ࠷খʹ͢ΔΑ͏ͳ#(supp(Q)) ≤ K Λຬͨ͢཭ࢄଌ౓Λ‫ٻ‬ΊΔ໰୊ʹରԠ͍ͯ͠Δ. ͜ͷ͜ͱʹ஫໨͢Δͱ, ୅ද఺ʹద੾ͳॏΈΛ༩͑Δ͜ͱͰ, KASP ͷ໰୊఺͕ղফ͞ΕΔ͜ͱ͕Θ͔Δ. ҰํͰ, ୅ද఺ͷ‫ݧܦ‬෼෍ͷζϨ͔Β, ฼ूஂ෼෍ʹ͓͍ͯີ౓ͷ௿͍఺΋୅ද఺ͱͯ͠ੜ੒͞ΕΔ. ͦΕΒͷ୅ද఺ʹ͸, ௿͍ॏΈׂ͕ΓৼΒΕΔͨΊ, ޮ཰͕ѱ͍.
͜ͷ໰୊Λղܾ͢ΔͨΊʹ, ൃදऀΒͷ͜Ε·Ͱͷ‫Ͱڀݚ‬, ৽͍͠୅ද఺ͷੜ੒ํ๏Ͱ͋Δ
Density-Preserving Vector Quantization (DPVQ) ΛఏҊ͍ͯ͠Δ. DPVQ ͕ੜ੒͢Δ୅ද఺ͷ‫ݧܦ‬෼෍͸, ઴ۙతʹσʔλͷഎ‫ޙ‬ͷ෼෍΁ऩଋ͢Δ͜ͱ͕ࣔͤΔ. ͦͷͨΊ, DPVQ ʹΑͬͯੜ੒͞Εͨ୅ද఺ʹ͸ฏ౳ͳॏΈׂ͕ΓৼΒΕΔͨΊ, ޮ཰తͳۙࣅ͕‫ظ‬଴Ͱ͖Δ. ҰํͰ, DPVQ͸, ີ౓ਪఆΛඞཁͱ͢ΔͨΊ, ࣍‫ͳ͘ߴ͕ݩ‬Δͱෆ҆ఆʹͳΔͱ͍͏໰୊͕͋Δ.

22

Algorithm 1 VQn (µ | r, K) ͷ࠷దԽΞϧΰϦζϜ
(0)

1: t ← 0 ͱ͠, Ϋϥελத৺ µ1(0), . . . , µK Λॳ‫ظ‬Խ͢Δɽ

2: for t = 0, . . . , T do
3:

!
(t) !
֤ i (i = 1, . . . , n) ʹରͯ͠ɼ!xi − µk ! Λ࠷খʹ͢ΔΫϥελ k ΛׂΓ౰ͯ, ‫ؼ‬ଐߦྻ"#(t)
U (t) = uij
ΛಘΔ.
n×K

4:


! !
!
1 if ∀j; !
!xi − µ(t) ! ≤ !xi − µ(t) !,
(t)
j
k
uik =
0 otherwise.

Ϋϥελฏ‫ۉ‬ΛҎԼͰߋ৽͢Δ.

1

n
(

(t+1)
(t) (t)µ̂k='uik wik xi ,(t) (t)n
j=1 ujk wjk i=1

!
!x − µ̂(t) !
!r−2

(t)
wik =


i

k

δ

͜͜Ͱ, δ > 0 ͸খ͍͞ਖ਼ͷఆ਺Ͱ͋Δ.
5:ऩଋ൑ఆ৚݅Λຬͨͤ͹ఀࢭ͠, ຬͨ͞ͳ͚Ε͹ t ← t + 1 ͱ͢Δ.
6: end for!
(t) !
if !xi − µ̂k ! > 0,!
(t) !
if !xi − µ̂ ! = 0.
kͦ͜Ͱ, ຊൃදͰ͸, ҎԼͰఆٛ͞ΕΔ order r ͷϕΫτϧྔࢠԽ‫ث‬Λ༻͍ͨۙࣅ๏Λߟ͑Δ.
n

VQn (µ | r, K) :=1(
min %xi − µk %r
1≤k≤K
n
i=1

͜͜Ͱ, r ∈ (0, 2] ͸ఆ਺, % · % ͸ Rd ্ͷϊϧϜ, x1 , . . . , xn ∈ Rd ͸֤σʔλ఺, µk ∈ Rd ͸ k ൪໨ͷΫϥελத৺, µ = (µ1 , . . . , µK ) Ͱ͋Δ. ຊൃදͰ͸, ‫ࢉܭ‬ͷ؆ศੑͷͨΊʹ, % · % ͸ Euclid ϊϧϜͱ͢Δ. r = 2 ͱ͢Ε͹ VQn (µ | r, K) ͸ k-means ๏ͱҰக͢Δ͕, r Λখ͘͞ͱΔ͜ͱͰ୅ද
఺ͷ‫ݧܦ‬෼෍ͱ฼ूஂ෼෍ͷζϨΛখ͘͢͞Δ͜ͱ͕Ͱ͖Δ. ຊൃදͰ͸, VQn (µ | r, K) ʹର͢Δ
࠷దԽ໰୊Λߴ଎ʹղͨ͘Ίʹ, k-means like ͳΞϧΰϦζϜΛఏҊ͢Δ. Algorithm 1 ʹ͓͍ͯ,VQn (µ(t) | r, K) > VQn (µ(t+1) | r, K) ͱ͍͏୯ௐ‫ݮ‬গੑ͕੒ΓཱͭͨΊ, ఀཹ఺΁ͷऩଋੑ͕อূ
Ͱ͖Δ.
ৄࡉͳཧ࿦తੑ࣭΍ r < 2 ͱͨ͠ϕΫτϧྔࢠԽ‫ث‬Λ༻͍ͨۙࣅͱ k-means ๏΍ DPVQ ʹΑΔ
ۙࣅͷ਺஋࣮‫ʹݧ‬ΑΔൺֱ͸౰೔ใࠂ͢Δ. ಛʹ, spectral clustering ʹରԠ͢Δ normalized cut(e.g., Terada and Yamamoto, 2019) ʹରͯ͠, ୅ද఺ͷ෼෍ͱ฼ूஂ෼෍ͷζϨͷӨ‫ڹ‬Λ໌Β͔ʹ͢Δ.
ࢀߟจ‫ݙ‬

Terada, Y. and Yamamoto, M. (2019). Kernel Normalized Cut: a Theoretical Revisit. In
Proceedings of the 36th International Conference on Machine Learning, PMLR 97, 6206–
6214.
von Luxburg, U. (2007). A tutorial on spectral clustering. Statistics and Computing, 17, 395–
416.
Yan, D., Huang, L., and Jordan, M. I. (2009) Fast approximate spectral clustering. In Proceedings of the 15th ACM SIGKDD International Conference on Knowledge Discovery and Data
Mining, 907–916.

23܏޲είΞղੳͷͨΊͷ৘ใྔ‫ن‬४͓Αͼࠩ෼ͷࠩ๏΁ͷల։ೋ‫ ٶ‬Յߦ
౷‫਺ܭ‬ཧ‫਺ ॴڀݚ‬ཧɾਪ࿦‫ܥڀݚ‬ຊߨԋͰ͸ɼ‫ج‬ຊతͳҼՌਪ࿦Ϟσϧ

y=

H
!
h=1

t(h) y (h) ,y (h) ∼ f (· | x(h) ; θ)Λѻ͏ɽ͜͜Ͱɼt(h) (∈ {0, 1}) ͸ h ൪໨ͷॲஔΛׂΓ౰ͯΒΕͨͱ͖ʹ 1 ͱͳΔׂΓ౰ͯม਺"(h) = 1ʣɼy (h) (∈ R) ͸ h ൪໨ͷॲஔΛׂΓ౰ͯΒΕͨͱ͖ͷજࡏ݁Ռม਺ɼx(h) (∈ Rr )
ʢ H
h=1 t
͸ y (h) ͷͨΊͷઆ໌ม਺ɼf (y (h) | x(h) ; θ) ͸ y (h) ʹର͢Δ x(h) ͷճ‫ؼ‬Ϟσϧɼθ (∈ Rp ) ͸ͦ͜Ͱ༻͍ΒΕΔύϥϝʔλͰ͋Δʢh ∈ {1, 2, . . . , H}ʣɽx(h) ͸ަབྷม਺ͷҰ෦Λ‫ؚ‬ΜͰ΋Α͍ɽ͜ͷϞσϧͰ͸ɼt(h) = 0 ͳΔ H − 1 ‫ݸ‬ͷ y (h) ͕ܽଌ͍ͯ͠Δͱߟ͑ΒΕΔɽͦͯ͠ɼҰൠʹ͸E[y (h) ] $= E[y (h) | t(h) = 1] Ͱ͋ΔͨΊɼ‫؍‬ଌ஋ͷΈ͔Β୯७ʹ θ Λਪఆ͢ΔͱόΠΞε͕ੜͯ͡͠·͏ɽͦ͜Ͱɼy (h) ͱ t(h) ͷަབྷม਺ z (∈ Rs ) ͕‫؍‬ଌ͞Ε͍ͯΔ͜ͱΛ૝ఆ͢Δɽͦͯ͠ɼ͜
ͷόΠΞεͷআ‫ڈ‬ΛՄೳͱ͢ΔͨΊͷ΋ͷͰ͋ΔɼແࢹͰ͖ΔׂΓ౰ͯ৚݅{y (1) , y (2) , . . . , y (H) } ⊥⊥{t(1) , t(2) , . . . , t(H) } | zΛԾఆ͢Δɽ·ͨɼਖ਼஋৚݅ P(t(h) = 1) > 0 ΛԾఆ͢Δɽ͜ͷϞσϧʹ͕ͨ͠͏ N ‫ݸ‬ͷαϯϓϧ
(1) (2)(H) (1)
͕͋Δͱ͠ɼୈ i αϯϓϧͷม਺ʹ͸ఴ͑ࣈ i Λ෇͚Δ͜ͱʹ͢Δɽui = (yi , yi , . . . , yi , ti ,
(2)
(H)
(1)
(2)
(H)
ti , . . . , ti , xi , xi , . . . , xi , zi ) ͱॻ͖ɼαϯϓϧ͸ಠཱɼͭ·Γ

ui ⊥
⊥ uj

(i $= j; i, j ∈ {1, 2, . . . , N })Ͱ͋Δ͜ͱ΋Ծఆ͢Δɽਪఆର৅͸ɼୈ k ॲஔ‫ࡍ࣮͕܈‬ΑΓ d(k) ഒଘࡏ͢Δ฼ूஂʹ͓͚Δύϥϝʔλ θ ͱ͢ΔɽҎ߱ɼύϥϝʔλͷਅ஋͸ θ ∗ ͷΑ͏ʹද͢͜ͱͱ͢Δɽy (h) ͱ z ͷؔ܎Λਖ਼͘͠ϞσϦϯάͰ͖Ε͹ɼ࠷໬๏͔ΒҼՌޮՌΛҰகਪఆ͢Δ͜ͱ͕Ͱ͖
Δ͕ɼҰൠʹ͜ͷϞσϦϯά͸೉͍͠ɽͦͷͨΊɼ͜ͷϞσϦϯάΛඞͣ͠΋ඞཁͱ͠ͳ͍ɼ܏޲
είΞ e(h) (z; α) ≡ P(t(h) = 1 | z; α) Λ༻͍ͨηϛύϥϝτϦοΫΞϓϩʔν͕༻͍ΒΕΔ͜ͱ
΋ଟ͍ɽ͜͜Ͱɼα (∈ Rq ) ͸܏޲είΞͷؔ਺Λಛ௃෇͚ΔύϥϝʔλͰ͋ΔɽͦͷΞϓϩʔν"(k) (k)(h) (z; α) ͱ͠ɼͷதͰͷ‫ج‬ຊ͸ɼ‫཰֬ٯ‬ॏΈ෇͚ਪఆͰ͋Δɽw(h) (z; α) ≡ H
k=1 d e (z; α)/e
N !
H
!
∂
(h)
(h)
(h)
ti w(h) (zi ; α)log f (yi | xi ; θ) = 0p∂θi=1 h=1Λղ͘͜ͱͰ‫཰֬ٯ‬ॏΈ෇͚ਪఆྔ θ̂ IPW ͸༩͑ΒΕΔɽແࢹͰ͖ΔׂΓ౰ͯ৚݅ͷ΋ͱɼ‫֬ٯ‬཰ॏΈ෇͚ਪఆྔ͸ҰகੑΛ΋ͭɽຊઃఆʹ͓͍ͯ y (h) ͸ z ͱ૬͕ؔ͋Δ͕ɼ‫཰֬ٯ‬ॏΈ෇͚ਪఆ͸ͦͷ‫ظ‬଴஋ͷਪఆʹ z ͷ৘ใΛ௚઀͸༻͍͍ͯͳ͍ɽೋॏ‫݈ؤ‬ਪఆ͸ɼͦΕΛ࣮ߦͯ͠‫཰֬ٯ‬ॏΈ෇͚ਪఆΛվྑͨ͠΋ͷͰ͋Δɽz Λ༩͑ͨ΋ͱͰͷ y (h) ͷ৚݅෇͖෼෍Λ͞ΒͳΔύϥϝʔλ β (∈ Rr ) Λ༻͍ͯp(h) (y (h) | z; β) ͱදͨ͠ͱ͖ɼlog f (y (h) | x(h) ; θ) Λ p(h) (y (h) | z; β) Ͱ‫ظ‬଴஋ΛͱΓɼͦͷ৚݅෇͖‫ظ‬଴஋Λ g (h) (x(h) , z; θ, β) ͱॻ͍ͯ༻͍Δɽ࣮ࡍʹ β ʹରͯ͠͸ɼͳΜΒ͔ͷҰகਪఆྔ β̂Λ୅ೖ͢Δɽ۩ମతʹ͸ɼೋॏ‫݈ؤ‬ਪఆྔ θ̂ DR ͸ɼ্‫ه‬ਪఆํఔࣜͷࠨลʹ
$
N !
H #!
H
!
∂ (h) (h)(h) (h)(k) (k)d ti − ti w (zi ; α)g (xi , zi ; θ, β)∂θ
i=1 h=1k=124ΛՃ͑ɼα ͱ β ʹਪఆྔΛ୅ೖ͠ɼͦΕΠίʔϧ 0p Λ θ ʹ͍ͭͯղ͘͜ͱͰ༩͑ΒΕΔɽ͜ͷਪఆํఔࣜͷ i ʹґଘ͢Δ෦෼Λ m(ui ; θ, α, β) ͱද͓ͯ͘͠ɽ܏޲είΞ͔৚݅෇͖‫ظ‬଴஋ͷͲͪΒ͔͑͞ਖ਼͘͠ϞσϦϯά͞Ε͍ͯΕ͹ θ̂ DR ͸ҰகੑΛ΋ͪɼͦΕ‫ʹނ‬ೋॏ‫ݺͱ݈ؤ‬͹ΕΔɽҼՌਪ࿦ʹ͓͚Δ৘ใྔ‫ن‬४Λಋͨ͘ΊͷϦεΫΛఆٛ͢Δɽ௨ৗͷ AIC λΠϓͷ৘ใྔ‫ن‬४(h)† (h)†(h)† †(h) (h)(h)Λߟ͑Δͱ͖ͱಉ༷ɼ(yi , ti , xi , zi ) Λ (yi , ti , xi , zi ) ͷίϐʔͱ͠ɼ‫཰֬ٯ‬ॏΈ෇͚ਪఆͰ༻͍ΒΕΔଛࣦؔ਺͔Βࣗવʹఆٛ͞ΕΔ

−2

N !
H
!

(h)†

E[ti

(h)†

w(h) (zi† ) log f (yii=1 h=1

(h)†

| xi

; θ̂ DR )]Λߟ͑Δɽ͜͜Ͱɼw(h) (z) ͸ w(h) (z; α̂) ͷऩଋઌͱ͢Δɽ͜ͷϦεΫ͸ɼୈ k ॲஔ‫ࡍ࣮͕܈‬ΑΓ
d(k) ഒଘࡏ͢Δ฼ूஂʹ͓͚Δɼਅͷ෼෍ͱਪఆ͞Εͨ෼෍ͷ Kullback-Leibler μΠόʔδΣϯ" "H (h) (h)(h)εʹ‫ͮ͘ج‬΋ͷͱΈͳͤΔɽ͜ͷϦεΫͷࣗવͳਪఆྔ͸ −2 N(zi ) log f (yi |i=1
h=1 ti w(h)xi ; θ̂ DR ) Ͱ͋Δ͕ɼͦΕ͸‫ݴ‬Θͣ΋͕ͳաখධՁ͢Δɽͦ͜ͰɼόΠΞεΛ
%!
&
N !
H
N !
H
!
(h) (h)(h)(h) DR(h)† (h) †(h)†(h)† DR−2E
ti w (zi ) log f (yi | xi ; θ̂ ) −ti w (zi ) log f (yi | xi ; θ̂ )i=1 h=1i=1 h=1ͱද‫͠ݱ‬ɼ͜ͷ‫ظ‬଴஋ͷத਎ͷऩଋઌΛ blimit ͱ͠ɼE(blimit ) Λ઴ۙόΠΞεͱͯ͠ิਖ਼ʹ༻͍
Δɽ઴ۙόΠΞεΛಋग़͢Δલʹɼ४උͱͯ͠ҎԼΛ༩͓͑ͯ͘ɽ

ิ୊. ʮ݁Ռม਺ʹର͢Δަབྷม਺ͷϞσϧʯ͔ʮׂΓ౰ͯม਺ʹର͢Δަབྷม਺ͷϞσϧʯͷͲͪΒ͔Ұํ͕ਖ਼͚͠Ε͹ɼೋॏ‫݈ؤ‬ਪఆྔͷ‫ࠩޡ‬͸
N #
H
!
1 !
(h) ∂
DR
∗
∗
† −1
∗
†
†
∗
†
†
θ̂ − θ =A(θ , α ) m(ui ; θ , α , β ) + C1 (θ , α , β )tilog e(h) (zi ; α† )
N
∂α
i=1
h=1
$
H
! ∂ (h)∗††(h) (h)†t log p (yi | zi ; β ) {1 + oP (1)}
+ C2 (θ , α , β )
∂β i
h=1

ͷ‫Ͱܗ‬ද͞ΕΔɽα† ɼβ † ͷͲͪΒ͔Ұํ͸ਅ஋ α∗ ɼβ ∗ Ͱ͋ΔɽAɼC1 ɼC2 ͷఆٛ͸͜͜Ͱ͸লུ͢Δ͕ɼC1 (θ ∗ , α† , β ∗ ) = C2 (θ ∗ , α∗ , β † ) = O ͕੒ཱ͍ͯ͠
Δɽͭ·Γɼ্‫ʹࠩޡه‬͸ଥ౰Ͱͳ͍Մೳੑͷ͋ΔείΞؔ਺͕‫·ؚ‬Ε͍ͯΔ͕ɼଥ౰Ͱͳ͍ͱ
͖͸܎਺͕ O ͱͳ͓ͬͯΓɼ໰୊͸ଘࡏ͠ͳ͍ɽ͜ͷิ୊ΑΓɼҎԼͷओ݁Ռ͕ಘΒΕΔɽ
ఆཧ. ʮ݁Ռม਺ʹର͢Δަབྷม਺ͷϞσϧʯ͔ʮׂΓ౰ͯม਺ʹର͢Δަབྷม਺ͷϞσϧʯͷͲͪΒ͔Ұํ͕ਖ਼͚͠Ε͹ɼೋॏ‫݈ؤ‬ਪఆͷͨΊͷ৘ใྔ‫ن‬४ͷ઴ۙόΠΞε͸

E(blimit ) = tr[A(θ ∗ , α∗ )−1 {B(θ ∗ , α∗ ) + D1 (θ ∗ , α∗ , β ∗ , α† , β † )} + D2 (θ ∗ , α∗ , β ∗ , α† , β † )]ͷ‫Ͱܗ‬༩͑ΒΕΔɽ

AɼB ɼD1 ɼD2 ͷఆٛ͸͜͜Ͱ͸লུ͢Δ͕ɼ‫ط‬ଘͷ৘ใྔ‫ن‬४͔Β૝૾Ͱ͖ΔΑ͏ʹɼ΍͸Γ઴ۙόΠΞε͸ύϥϝʔλͷਅ஋ʹґଘ͍ͯ͠Δ͜ͱ͕Θ͔Δɽα∗ ͔ β ∗ ͷͲͪΒ͔ΛਪఆͰ͖ͳ͍ຊઃఆͰ͸ɼ઴ۙόΠΞεͷධՁ͸Ͱ͖ͳ͍ͱࢥΘΕΔ͔΋͠Εͳ͍͕ɼਅ஋Λਪఆͤͣʹ‫ܦ‬‫ݧ‬ਪఆΛ͏·͘༻͍Δ͜ͱͰɼD1 ͱ D2 ͷҰகਪఆྔ D̂1 ͱ D̂2 Λ༩͑Δ͜ͱ͕Ͱ͖Δɽ݁ՌɼDRIC ≡ − 2

N !
H
!

(h)

ti
(h)
(h)
log f (yi | xi ; θ̂ DR )
(h)
e
(z
;
α̂)
i
i=1 h=1+ 2tr[A(θ̂ DR , α̂)−1 {B(θ̂ DR , α̂) + D̂1 (θ̂ DR , α̂, β̂)} + D̂2 (θ̂ DR , α̂, β̂)]ͱ͍͏‫ܗ‬ͷ৘ใྔ‫ن‬४Λ༩͑Δ͜ͱ͕Ͱ͖Δɽ৘ใྔ‫ن‬४ࣗମ͕ೋॏ‫݈ؤ‬ੑΛ΋͍ͬͯΔͨΊɼೋॏ‫ج݈ؤ‬४ͱ‫͢ʹͱ͜Ϳݺ‬Δɽ

25

外れ値に頑健な因果効果の推定量
原 田 和 治 a ，藤 澤 洋 徳 ba) 東 京 医 科 大 学 医 療 デ ー タ サ イ エ ン ス 分 野b) 統 計 数 理 研 究 所 数 理・推 論 研 究 系

1 はじめに

統 計 学 的 推 定 に お い て ，外 れ 値 は バ イ ア ス の 主 な 原 因 の ひ と つ で あ る ．Canavire-Bacarrez

ら は 平 均 因 果 効 果 の 様 々 な 推 定 量 に つ い て ，外 れ 値 の 影 響 を 実 験 的 に 評 価 し た [1]．特 に ，vertical outlier と 呼 ば れ る ，結 果 変 数 に 入 る タ イ プ の 外 れ 値 が ，平 均 因 果 効 果 の 推 定 に 深 刻なバイアスをもたらすことが指摘されている．

ロ バ ス ト 統 計 は ，外 れ 値 に 頑 健 な 統 計 学 的 推 論 の 方 法 を 扱 う 分 野 で あ る ．平 均 因 果 効 果に 対 す る ロ バ ス ト な 推 定 量 と し て は ，例 え ば ，中 央 値 因 果 効 果 に 対 す る 推 定 量 が 挙 げ ら れる ．具 体 的 に は ，Zhang ら の IPW お よ び 二 重 頑 健（Doubly Robust; DR）中 央 値 な ど が あ る[2]．平 均 因 果 効 果 と は 推 定 対 象 が 異 な る も の の ，後 述 す る 潜 在 ア ウ ト カ ム が 対 称 分 布 に 従
う な ら ば ，中 央 値 因 果 効 果 に 対 す る 推 定 量 を 平 均 因 果 効 果 に 対 す る ロ バ ス ト 推 定 量 と して 用 い る こ と が で き る ．し か し ，一 般 に 中 央 値 は 平 均 に 比 べ て 外 れ 値 に 頑 健 で は あ る も のの ，影 響 を 取 り 除 け る わ け で は な い ．特 に ，外 れ 値 が 分 布 の 中 心 か ら 見 て 片 方 に 偏 っ て い
る 場 合 に ，外 れ 値 の 割 合 が 高 ま る ほ ど ，そ の 影 響 は 顕 著 と な る ．

本 発 表 で は ，DR 推 定 量 を 拡 張 し ，平 均 因 果 効 果 に 対 す る 中 央 値 よ り も 外 れ 値 に 頑 健 な

推 定 量 を 提 案 し た ．提 案 手 法 の 基 本 的 な ア イ デ ア は ，中 央 値 よ り も 頑 健 な 密 度 べ き 重 み づ

け 推 定 方 程 式 に ，二 重 頑 健 な M 推 定 量（DR M 推 定 量 ）の 枠 組 み を 適 用 す る も の で あ る ．
し か し ，単 に DR M 推 定 量 を 当 て は め る だ け で は ，モ デ ル 選 択 に 関 す る 二 重 頑 健 性 を 失 う
こ と を 示 す ．さ ら に ，外 れ 値 の 割 合 に よ る 補 正 に よ り ，二 重 頑 健 性 を 保 ち な が ら ，外 れ 値に も 頑 健 な 推 定 量 を 構 成 で き る こ と を 示 し た ．な お ，本 発 表 は 査 読 付 き 学 術 誌 に 掲 載 さ れた 論 文 [3] を も と に し て い る ．

2 提案手法提 案 し た ３ 手 法（DP-IPW, DP-DR, εDP-DR）の う ち ，二 重 頑 健 性 を 保 ち な が ら ，外 れ 値に 対 し て も 頑 健 な ε DP-DR 法 の 概 略 を 示 す ．

DP-DR 法 は ，DR 推 定 量 を 一 般 の M 推 定 に 拡 張 し た DR M 推 定 量 に お い て ，推 定 関 数 に

密 度 べ き 型 推 定 関 数 を 選 択 し た も の で あ っ た ．し か し ，DP-DR 法 は 外 れ 値 存 在 下 で 二 重 頑
健 性 を 保 て な い こ と を 報 告 し た ．DP-DR 法 の 外 れ 値 存 在 下 で の バ イ ア ス を 補 正 し ，外 れ 値

存 在 下 で も 二 重 頑 健 性 を 満 た す よ う に 構 成 し た も の が ，εDP-DR 推 定 量 で あ る ．次 の 推 定

26

方 程 式 の 解 を εDP-DR 推 定 量 と 定 義 す る ．
n

1!
n i=1

"

Ti
Ti − π(Xi ; α̂)h(Yi ; µ)γ (Yi − µ) − (1 − ε̂)Eq̂ [h(Yi ; µ)γ (Yi − µ)|Ti = 1, Xi ]π(Xi ; α̂)π(Xi ; α̂)

#

=0
(1)

外 れ 値 の 割 合 ε は ，例 え ば 先 述 の [4] に よ っ て ，条 件 付 分 布 q と 同 時 に 一 致 推 定 で き る ．

εDP-DR 推 定 量 は ，外 れ 値 が な い 場 合 に は DP-DR 推 定 量 と 同 様 に 二 重 頑 健 性 を 持 つ ．さ らに ，外 れ 値 が あ る 場 合 で も ，推 定 方 程 式 の 不 偏 性 の 意 味 で 二 重 頑 健 性 を 保 持 し て い る ．傾向 ス コ ア が 正 し い 場 合 の 証 明 は DP-DR 推 定 量 と 同 様 で あ り ，条 件 付 分 布 モ デ ル が 正 し い

場 合 に も ，補 正 に よ り ，推 定 方 程 式 の 不 偏 性 を 証 明 で き る ．ε̂ に よ る 補 正 は 傾 向 ス コ ア が正 し い 場 合 は 必 要 な く ，条 件 付 分 布 モ デ ル だ け が 正 し い 場 合 に は ，ε̂ も 同 時 に 正 し く 推 定
さ れ る た め ，補 正 の た め に ，新 た な モ デ ル を 特 定 す る 必 要 は な い ．

3 数値実験

提 案 手 法 の 外 れ 値 耐 性 や 推 定 効 率 に つ い て ，実 験 的 に 検 討 し た ．性 能 評 価 の た め ，外 れ

値 耐 性 な し の IPW/DR 推 定 量 の ほ か ，各 種 提 案 さ れ て い る IPW/DR 中 央 値 と 比 較 し た ．
本 発 表 の 設 定 に 沿 っ た シ ミ ュ レ ー シ ョ ン 用 デ ー タ セ ッ ト を ，サ ン プ ル サ イ ズ n = 100 で10,000 セ ッ ト 生 成 し た ．真 の µ(1) は 3 と し た ．

全 般 的 な 傾 向 と し て ，提 案 手 法 は 外 れ 値 耐 性・推 定 効 率 の 両 面 で 中 央 値 よ り も 優 れ る 傾

向 に あ っ た ．ま た ，γ の 値 を 大 き く す る と 外 れ 値 耐 性 は 増 す も の の ，推 定 効 率 は 低 下 す る
傾 向 に あ っ た ．こ の 傾 向 は 他 の 密 度 べ き 重 み を 用 い た ロ バ ス ト 推 定 量 と 同 様 で あ る ．

ま た ，本 稿 が 想 定 す る 混 合 分 布 型 の 外 れ 値 と は 別 に ，誤 差 項 を Cauthy 分 布 と す る こ とで ，対 称 に 外 れ 値 が 発 生 す る 場 合 に つ い て も 検 討 し た ．外 れ 値 が 対 称 の 場 合 ，中 央 値 も 一
致 性 を 失 わ な い た め ，上 述 の 結 果 に 比 べ る と 中 央 値 に 有 利 な 設 定 で あ る ．そ れ に も 関 わ らず ，提 案 手 法 は 中 央 値 に 基 づ く 各 手 法 よ り も RMSE が 小 さ か っ た ．

加 え て ，DP-DR 法 と εDP-DR 法 を 比 較 し ，ε に よ る 補 正 の 効 果 を 検 証 し た 結 果 に つ い ても 報 告 し た ．条 件 付 分 布 モ デ ル だ け を 正 し く 特 定 し ，各 手 法 を 適 用 し た 結 果 を 表 3 に 示 す ．

DP-DR 推 定 量 は ，γ を 十 分 大 き く し て も ，外 れ 値 の 増 加 に 伴 い バ イ ア ス が 拡 大 し て い る ．
一 方 ，εDP-DR 法 は つ ね に 真 値 に 近 い 値 を 推 定 で き て い る ．こ の 結 果 か ら も ，ε̂ に よ る 補 正がバイアスの軽減に機能していると考えられる．

参考文献[1] Gustavo Canavire-Bacarreza, Luis Castro Peñarrieta, and Darwin Ugarte Ontiveros. Outliers in
Semi-Parametric estimation of treatment effects. Econometrics, 9(2):19, April 2021.
[2] Zhiwei Zhang, Zhen Chen, James F Troendle, and Jun Zhang. Causal inference on quantiles withan obstetric application. Biometrics, 68(3):697–706, September 2012.
[3] Kazuharu Harada and Hironori Fujisawa. Outlier-Resistant estimators for average treatment effectin causal inference. Statistica Sinica, 34(2), April 2024. (to appear).
[4] Takafumi Kanamori and Hironori Fujisawa. Robust estimation under heavy contamination using
unnormalized models. Biometrika, 102(3):559–572, May 2015.

27

Statistical Inference for Local Granger Causality
Yan Liu, Masanobu Taniguchi, Hernando Ombao 1
Abstract
Granger causality has been employed to investigate causality relations between components ofstationary multiple time series. We generalize this concept by developing statistical inference
for local Granger causality for multivariate locally stationary processes. Our proposed local
Granger causality approach captures time-evolving causality relationships in nonstationary processes. The proposed local Granger causality is well represented in the frequency domain andestimated based on the parametric time-varying spectral density matrix using the local Whittle
likelihood. Under regularity conditions, we demonstrate that the estimators converge to multivariate normal in distribution. Additionally, the test statistic for the local Granger causality isshown to be asymptotically distributed as a quadratic form of a multivariate normal distribution. For practical demonstration, the proposed local Granger causality method uncovered newfunctional connectivity relationships between channels in brain signals. Moreover, the method
was able to identify structural changes in financial data.
keywords: Brain signals, Local Granger causality, Local Whittle likelihood, Multivariate locallystationary processes, Time-varying spectral density matrix, Topological data analysis
1. Local Granger Causality
[1]

[p]

Let Xt,T = (Xt,T , . . . , Xt,T )! be a sequence of p-dimensional multivariate stochastic processes
Xt,T =∞!

At,T (j)!t−j ,

(1.1)

j=−∞

where the sequences {At,T (j)}j∈Z satisfy the regularity conditions. The process (1.1) is usuallyreferred to as the multivariate locally stationary process.
" (1) !
(2) ! #!
Let m and M be two positive integers such that p = m+M . Suppose Xt,T = Xt,T , Xt,T
,
(1)

(2)

Xt,T ∈ Rm , Xt,T ∈ RM , has the time-varying spectral density matrix f (u, λ) with the partition$%1
f (u, λ)11 f (u, λ)12f (u, λ) =:=A(u, λ)KA(u, −λ)! ,u ∈ [0, 1],f (u, λ)21 f (u, λ)22
2π&∞where A(u, λ) :=j =−∞ A(u, j) exp( ijλ). Let Σ(u) be the one-step-ahead prediction error
covariance matrix based on the time-varying spectral density matrix f (u, λ) with the same
partition.
1 this work was supported by JSPS Grant-in-Aid for Scientific Research (C) 20K11719 (Liu, Y.) and JSPS

Grant-in-Aid for Scientific Research (S) 18H05290 (Taniguchi, M.)

28

Y. Liu et al. / Local Granger Causality
(1)

2

(2)

Let H(τ1 , τ2 ) = sp(Xt,T , 1 ≤ t ≤ τ1 ; Xt,T , 1 ≤ t ≤ τ2 ) be the closed linear subspace generated
(1)

(2)

by {Xt,T , 1 ≤ t ≤ τ1 ; Xt,T , 1 ≤ t ≤ τ2 }. Especially, we use H(τ1 , 0) and H(0, τ2 ) to express the
(1)

(2)

closed linear subspace generated by {Xt,T , t ≤ τ1 }, {Xt,T , t ≤ τ2 }, respectively. Introducing thecompanion process
" (2)
#
(2)
(2)
Yt,T = Xt,T − E Xt,T | H(t, t − 1) ,
(2)

(1)

we propose the local Granger causality measure from {Xt,T } to {Xt,T } as' π
1
GC(2→1) (u) =FGC(u, λ)dλ,
2π −π
where

|f (u, λ)11 |(.
FGC(u, λ) = ln ((((f (u, λ)11 − 2πg(u, λ)12 Σ̃(u)−122 g(u, λ)21 (" (1) ! (2) ! #!
Here, g(u, λ) is the time-varying spectral density matrix of the process { Xt,T , Yt,T}, andΣ̃(u) is an (M × M )-matrixΣ̃(u)22 = Σ(u)22 − Σ(u)21 Σ(u)−111 Σ(u)12 .
(2→1)

Regarding the local hypothesis H0: GC(2→1) (u) = 0, we use the following test statistic:S † (u) := 2T bT GC(2→1) (u; θ̂T ),where θ̂T is the Whittle estimate. Then we have the following result.
(2→1)

Theorem 1.1. Under the null" hypothesis# H0: GC(2→1) (u) = 0, if assumptions in Liu−1−6
et. al. (2021) hold with bT = o T (ln T )and bT = o(T −1/5 ), we have
"
#!
"
#
L
S † (u) −→ N 0, V(u) H(u)N 0, V(u) ,"#where N 0, V(u) and H(u) are defined in Liu et. al. (2021).

References

Dahlhaus, R. (2000). A likelihood approximation for locally stationary processes. The Annalsof Statistics 28 1762–1794.
Dahlhaus, R. (2009). Local inference for locally stationary time series based on the empirical
spectral measure. Journal of Econometrics 151 101–112.
Granger, C. W. J. (1963). Economic processes involving feedback. Information and Control
6 28–48.
Granger, C. W. J. (1969). Investigating causal relations by econometric models and crossspectral methods. Econometrica: Journal of the Econometric Society 37 424–438.
Hosoya, Y. (1991). The decomposition and measurement of the interdependency between
second-order stationary processes. Probability Theory and Related Fields 88 429–444.
Liu, Y., Taniguchi, M. and Ombao, H. (2021). Statistical inference for local Granger causality. arXiv: https://arxiv.org/abs/2103.00209.
Ombao, H., von Sachs, R. and Guo, W. (2005). SLEX analysis of multivariate nonstationary
time series. Journal of the American Statistical Association 100 519–531.
Taniguchi, M., Puri, M. L. and Kondo, M. (1996). Nonparametric approach for nonGaussian vector stationary processes. Journal of Multivariate Analysis 56 259–283.

29ऑ͍‫ڞ‬௨αϙʔτ৚݅ԼͰͷ֬཰తհೖʹ‫ͮ͘ج‬ҼՌޮՌͷਪఆࢁຊ ྙੜ 1,2
1 େࡕେֶେֶӃਓؒՊֶ‫ڀݚ‬Պ, 2 ཧԽֶ‫৽ֵॴڀݚ‬஌ೳ౷߹‫ڀݚ‬ηϯλʔ

1. ͸͡Ίʹ ࣾձ੓ࡦ΍‫ݸ‬ਓͷߦಈʹؔ͢ΔհೖͳͲɼܾఆ࿦తͳॲஔͷׂ෇͕‫࣮ݱ‬తͰͳ͍৔
߹ʹରԠ͢ΔͨΊͷํ๏ͱͯ֬͠཰తհೖʢstochastic interventionsʣ͕͋Γɼۙ೥ଟ͘ͷ‫͕ڀݚ‬ߦΘΕ͍ͯΔɽ֬཰తհೖͰ͸ɼԾ૝తͳॲஔͷׂ෇͸‫ڞ‬มྔʹґଘͯ֬͠཰తʹܾఆ͞ΕΔͱ
Ծఆ͢Δɽྫ͑͹ɼDı́az and van der Laan (2012) Ͱ͸ɼॲஔͷׂ෇ϝΧχζϜͱͯ͠հೖ෼෍
Λ૝ఆ͠ɼॲஔͷਫ४Λγϑτͤͯ֬͞཰తʹॲஔΛׂΓ෇͚Δํ๏ΛఏҊ͓ͯ͠Γɼ͜Ε·Ͱʹ͍͔ͭ͘ͷ֦ு͕ͳ͞Ε͍ͯΔɽҰൠతͳҼՌਪ࿦ख๏ͱಉ༷ʹɼDı́az and van der Laan (2012)ͷํ๏Ͱ΋ޮՌͷࣝผ΍ਪఆͷͨΊʹ͸ਖ਼஋ੑͷԾఆ͕ॏཁͱͳΔɽ͔͠͠ɼ֬཰తհೖͰѻ͏ଟ஋΍࿈ଓతͳਫ४Λ࣋ͭॲஔͰ͸ͦͷ੒ཱΛԾఆ͢Δͷ͸ඇ‫࣮ݱ‬తͰ͋Δɽͦ͜Ͱຊ‫Ͱڀݚ‬͸Sufficient Dimension Reduction (SDRʣʹ‫͍ͯͮج‬ਖ਼஋ੑʹର͢Δ৚݅Λ؇࿨ͨ͠ਪఆྔΛఏҊ͢Δɽ
2. ఏҊํ๏ X ∈ X ⊂ Rp Λ‫ڞ‬มྔϕΫτϧɼT ∈ T Λॲஔม਺ɼY ∈ R ΛΞ΢τΧϜͱ͠ɼ(X, T, Y ) ͸෼෍ P0 ʹै͏ͱ͢Δɽ·ͨɼಠཱಉ෼෍ͳ‫؍‬ଌσʔλ (X i , Ti , Yi ) ∼ P0 (i = 1, . . . , n)͕ಘΒΕ͍ͯΔͱ͢ΔɽϊϯύϥϝτϦοΫߏ଄ํఔࣜϞσϧ X = fX (EX ), T = fT (X, ET ), Y =fY (T, X, EY ) ʹ͓͍ͯɼॲஔม਺ T ͸հೖ෼෍ Pδ (g0 )(T = t|X) = g0 (t − δ(X)|X) ʹैͬͯ஋͕ܾఆ͞ΕΔͱԾఆ͢Δɽ͜͜Ͱɼδ(X) ͸ॲஔͷਫ४Λγϑτͤ͞ΔྔΛද͠ɼg0 ͸ਅͷʢॲஔ΁ͷʣമ࿐ϝΧχζϜ g0 (t|X) = P0 (T = t|X) Λද͢ɽհೖ෼෍Λಋೖͨ͠΋ͱͰͷજࡏΞ΢τΧϜ Y (Pδ ) ͷ‫ظ‬଴஋ E[Y (Pδ )] ͕‫ڵ‬ຯͷ͋ΔύϥϝʔλͰ͋Δɽͳ͓ɼհೖʹΑΓ֬཰ 1 Ͱ T = tʹ‫ݻ‬ఆ͞Εͨ৔߹ͷજࡏΞ΢τΧϜ͕ɼҰൠతͳજࡏΞ΢τΧϜ Y (t) ʹରԠ͍ͯ͠Δɽ৚݅෇͖
ަ‫׵‬Մೳੑ Y (t)⊥⊥T | X, ∀t ∈ T ͓ΑͼҰகੑͷԾఆ T = t ⇒ Y = Y (t) ͷ΋ͱͰɼྫ͑͹ X ͕཭ࢄͷ৔߹ɼ‫ڵ‬ຯͷ͋Δύϥϝʔλ͸֬཰෼෍ P ͷؔ਺ Ψ(P ) ʹΑΓ! !
E[Y (Pδ )] = Ψ(P ) =E(Y |T, X)Pδ (g)(T |X)P (X)(1)
T ∈T X∈XͷΑ͏ʹද͞ΕΔɽͳ͓ɼύϥϝʔλͷਅ஋͸ਅͷ෼෍ P0 Λ༻͍ͯ ψ0 = Ψ(P0 ) ͱද͞ΕΔɽύϥϝʔλ (1) ͷਪఆʹ͓͍ͯɼDı́az and van der Laan (2012) Ͱ͸‫཰֬ٯ‬ॏΈ෇͚ (IPTW)ਪఆྔɼ֦ு IPTW (AIPTW) ਪఆྔɼ͓Αͼ targeted maximum likelihood estimator (TMLE)͕ఏҊ͞Ε͍ͯΔɽͦΕΒͷਪఆྔͷଥ౰ੑΛอূ͢ΔͨΊʹ͸ਖ਼஋ੑ৚͕݅ඞཁͰ͋ΓɼͦͷͨΊʹ͸‫ڞ͍ڧ‬௨αϙʔτ৚݅ʢLuo et al., 2017ʣ

Ω(X) = Ω(X | T = t), t ∈ T͕ඞཁͱͳΔɽͨͩ͠ɼΩ(X) ͓Αͼ Ω(X | T = t) ͸ͦΕͧΕ฼ूஂશମͰͷ X ͷऔΓ͏Δ஋ɼ
T = t ͷूஂͰͷ X ͷऔΓ͏Δ஋Λද͢ɽ͜ͷ৚݅͸֤ॲஔ‫Ͱ܈‬ͷ‫ڞ‬มྔ X ͷ෼෍͕‫׬‬શʹॏͳ͍ͬͯΔ͜ͱΛཁ‫͍ͯ͠ٻ‬Δ͕ɼX ͕ߴ࣍‫ݩ‬ͷ৔߹΍‫࣌ܦ‬ଌఆσʔλͰ࣌఺͝ͱʹ͜ͷ৚͕݅ඞཁͱ͞ΕΔ৔߹ͳͲɼͦͷ੒ཱΛԾఆ͢Δͷ͕೉͍͠৔߹͕ଟ͍ɽͦ͜ͰɼఏҊख๏Ͱ͸ SDR
Λ༻͍ͯΑΓऑ͍৚݅ͷ΋ͱͰύϥϝʔλ (1) ͷਪఆΛߦ͏ɽఏҊख๏Ͱ͸ɼ·֤ͣ t ∈ T Ͱͷજࡏ݁Ռม਺ʹର͢Δ Central Mean Subspace (CMS; Cookand Li, 2002) SE(Y (t)|X) Λ‫؍‬ଌσʔλ͔Βਪఆ͢Δඞཁ͕͋ΔɽҰൠʹ͸ɼ‫؍‬ଌσʔλ͔Βਪఆ30DՄೳͳ‫؍‬ଌΞ΢τΧϜ Y ʹ͓͚Δॲஔͷ֤ਫ४ t ͷ෦෼ूஂͰͷ CMS SE(Y|X,T =t) ͱ্‫ه‬ͷ CMS
SE(Y (t)|X) ͕ಉ͡Ͱ͋Δอূ͸ͳ͍ɽ͜ͷ໰୊ʹରͯ͠ɼఏҊख๏ͷઃఆͰ͸ɼLuo et al. (2017)ͷ Theorem 1 ͱಉ༷ʹɼ৚݅෇͖ަ‫׵‬ՄೳੑͱҰൠͷ T ʹର͢Δऑ͍‫ڞ‬௨αϙʔτ৚݅"Ω(B"t X) = Ω(Bt X | T = t), t ∈ Tͷ΋ͱͰɼೋͭͷ CMS ͕Ұக͢Δ͜ͱ͕ࣔ͞ΕΔɽΑͬͯɼ‫؍‬ଌσʔλ͔ΒਪఆՄೳͳ CMS Λɼજࡏ݁Ռม਺ʹର͢Δ CMS ͱಉҰࢹͯ͠ར༻͢Δ͜ͱ͕Ͱ͖ΔɽఏҊख๏Ͱ͸ɼ·֤ͣ t ∈ T ʹର͢Δ CMS SE(Y (t)|X) Λ minimum average variance estimator(Xia et al., 2002) Λ༻͍ͯਪఆ͠ɼͦͷ‫ج‬ఈϕΫτϧ͔ΒͳΔ܎਺ߦྻ Bt Λਪఆ͢Δɽ۩ମతʹ͸ɼLuo et al. (2017) ͱಉ༷ʹɼ֤ t ∈ T ʹରͯ͠ɼn !
n!
i=1 j=1"#
I(Ti = t){Yi − at (X j , Bt ) − (X i − X j )" Bt bt (X j )}2 Kh1 ,r(t) B"
t (X i − X j )(2)Λ࠷খʹ͢Δ Bt ∈ Rp×r(t) Λ‫ٻ‬ΊΔɽ͜͜ͰɼKh,r (·) = hr K(·/h) ͸όϯυ෯ h Λ࣋ͭ Rr ্ͷΧʔωϧؔ਺Ͱ͋Δɽ(2) Λ࠷খʹ͢Δ at (x, Bt )ɼbt (x)ɼBt ͸ަ‫࠷ޓ‬খ 2 ৐๏ʹΑͬͯ‫ٻ‬ΊΔ͜$ t ͱ͢Δɽ
ͱ͕Ͱ͖Δɽਪఆ͞Εͨ܎਺ߦྻΛ B$ " X) Λ࣍ʹɼ֤ t ∈ T ʹରͯ͠ɼΞ΢τΧϜͷ৚݅෇͖‫ظ‬଴஋ E(Y |T, B
T$ t)
Wh ( B

n
!
i=1

"
#2
"
#
$ t ) − (X i − x)" B$ t b$ (x) Kh ,r(t) B̂" (X i − x)I(Ti = t) Yi − a$t (x, B
t
t
2

(3)

Λ࠷খʹ͢Δ a$t (x, B̂t ) ʹΑͬͯਪఆ͢Δɽͨͩ͠ɼCMS ͷ܎਺ߦྻ Bt ʹରͯ͠ Wh (Bt ) =%−1 Ͱ͋Δɽͳ͓ɼ(3) Λ࠷খʹ͢Δ a$ (x, B̂) ͓Αͼ b$ (x) ∈ Rr(t) ͸ɼ[ ni=1 Kh,r(t) {B"t (X i − x)}]tt$ t ) ͷਪ௨ৗͷ‫ॴہ‬ઢ‫ܗ‬ճ‫ͱؼ‬ಉ༷ʹॏΈ෇͖࠷খ̎৐ਪఆྔͱͯ͠‫ٻ‬ΊΔ͜ͱ͕Ͱ͖Δɽa$t (x, B$ t ) ͱ͢ΔɽఆྔΛ â$t (x, B࠷‫ʹޙ‬ɼ‫ڵ‬ຯͷ͋Δύϥϝʔλ (1) Λn1!
$ T +δ(X ) )I(Ti + δ(X i ) ∈ T )âTi +δ(X i ) (X i , B
i
i
n

(4)

i=1

ʹΑͬͯਪఆ͢Δɽ
CMS ͷਪఆ͓ΑͼΞ΢τΧϜͷ৚݅෇͖‫ظ‬଴஋ͷਪఆʹؔ͢Δ͍͔ͭ͘ͷԾఆͷ΋ͱͰਪఆྔ
(4) ͸Ұகੑ͓Αͼ઴ۙਖ਼‫ن‬ੑΛ࣋ͭ͜ͱ͕ࣔ͞ΕΔɽγϛϡϨʔγϣϯʹΑΓɼಛʹ‫ڞ͍ڧ‬௨αϙʔτ৚͕݅੒ཱͤͣɼΞ΢τΧϜͱ‫ڞ‬มྔɼ͓Αͼɼॲஔͱ‫ڞ‬มྔͷؔ܎͕ඇઢ‫͋Ͱܗ‬Δ৔߹ʹɼ‫ط‬ଘͷख๏ΑΓ΋ఏҊख๏͕ྑ͍ੑ࣭Λ࣋ͭ͜ͱ͕֬ೝ͞Εͨɽࢀߟจ‫ݙ‬

Cook, R. D. and Li, B. (2002). Dimension reduction for conditional mean in regression. Annalsof Statistics, 30, 455–474.
Dı́az, I., and van der Laan, M. J. (2012). Population intervention causal effects based onstochastic interventions. Biometrics, 68, 541ʕ549.
Luo, W., Zhu, Y., and Ghosh, D. (2017). On estimating regression-based causal effects using
sufficient dimension reduction. Biometrika, 104, 51–65.
Xia, Y., Tong, H., Li, W. K., and Zhu, L.-X. (2002). An adaptive estimation of dimension
reduction space (with discussion). Journal of the Royal Statistical Society, Series B, 64,
364–410.

31