f:id:kur0cky:20201023211747j:plain

こんにちは，DSOC 研究員の黒木裕鷹です．今年の4月に新卒として入社しましたが，あっという間に半年経ってしまいました．新しい環境にも慣れてきたのですが，同時に時間の速さに焦っています．

この連載では，私の学習や復習もかねて，ネットワーク構造を持ったデータ（グラフデータ）を対象とした統計解析を紹介しています．第2回のテーマも前回に引き続き，ネットワークデータの標本調査です．第1回で扱った，単純なノードやエッジの無作為抽出に加えて，第2回では少し凝った標本抽出を2つ紹介します．また，それらの標本抽出を使った，観測しきれないネットワーク全体の大きさ（ノードの数）を推定する取り組みを紹介します．それでは，やっていきます．

前回のおさらい
その他の標本抽出
- snowball sampling
- link tracing
ネットワークの大きさを推定する
- link tracing を用いた物理インターネットの大きさの推定 (Viger et al., 2007)
- snowball sampling に基づいた，ネットワークの大きさの推定 (Frank and Snijders, 1994)
おわりに
参考文献

前回のおさらい

buildersbox.corp-sansan.com

前回記事では，標本調査の初歩の初歩と，ネットワークの標本抽出を扱いました．ネットワークは構造をもつデータであるため，様々な標本抽出を考えることができます．中でも最も単純なものである induced sampling と incident sampling を紹介しました． induced sampling はノード単位で，incident sampling はエッジ単位で無作為抽出をそれぞれ行い，それらに基づいた部分グラフを抽出します．それぞれのイメージを説明していきます．

いま，図1のようなネットワーク構造からの標本抽出を考えます．ここでは説明のため小さいネットワークを作りましたが，人間の交友関係やWorld Wide Web (WWW) のように，巨大なネットワークは通常その全貌を観測することができません．そこで，標本抽出と統計的な推論によってその性質を測ることになります．

f:id:kur0cky:20201015100743p:plain — 図1：母集団ネットワーク

induced sampling のイメージを図2に示しました．まず，ノードの無作為抽出により，黄色のノードが観測されたとします．続いて，観測されたノード間のエッジ（太線）が観測されることになります．

f:id:kur0cky:20201015100746p:plain — 図2：induced sampling

incident sampling の例を図3に示しました．まず，エッジの無作為抽出により，太線のエッジが観測されたとします．続いて，観測されたエッジに紐づいているノード（黄色）が観測されることになります．

f:id:kur0cky:20201015100757p:plain — 図3：incident sampling

その他の標本抽出

前回紹介した2つの単純な標本抽出に加え，新たに2つのネットワークの抽出，snowball sampling と link tracing を紹介します．これらはいずれも，ネットワークの構造を積極的に利用した抽出で，実際の調査で用いられることも多いです． snowball sampling におけるノードの包含確率など，より詳細な性質についてはネットワーク解析の代表的な教科書である "Statistical Analysis of Network Data" (Kolaczyk, 2009)などを参照ください．

snowball sampling

snowball sampling の手順を箇条書きで簡単に説明します．また，そのイメージを図4に示しました．

seed とよばれる $N _ {(0)} ^ *$ 個の初期ノード（図では黄色）を無作為抽出し，その induced subgraph $\mathcal G ^ * _ {(0)}$ を観測します．
次に，seed の近傍のノード（オレンジ），seed との間のエッジ（太線）も観測し， $\mathcal G ^ * _ {(1)}$ を得ます．
この作業を $w$ 回繰り返し，得られた $\mathcal G ^ * _ {(w)}$ を観測ネットワーク $\mathcal G ^ *$ とします．

近傍，さらにその近傍，といったように，seed から雪玉のように大きく膨らませて抽出するため，snowball sampling という名前がついています．イメージでは $w = 1$ までを示していますが， $w = 2$ の snowball sampling を行えば，この例では母集団ネットワークの全貌を観測できることになります．

f:id:kur0cky:20201015104346p:plain — 図4：snowball sampling

snowball sampling は，社会調査で頻繁に用いられます．たとえば違法薬物の使用者やホームレスなどは，研究者にとってアクセスするのが難しくなるでしょう．ある専門知識を持った集団などを対象とする場合も，その存在は稀有であるため，無作為抽出では観測することが困難です．しかし，このような隠れた母集団の構成員は互いに密接に関わっていることが想定されるため，snowball sampling が役に立ちます．一方で，snowball sampling は seed の影響を強く受けてしまうことがデメリットとして挙げられます．これは，seed の属するコミュニティに偏った標本が観測されやすくなるためです．

link tracing

link tracing では，あるノード集合と別のノード集合を用意し，その経路上にあるノードとエッジを抽出します．通常，経路には最短経路（のうちの1つ）が使われることが多いようです．起点となるノード集合はソースノード，終点となるノード集合はターゲットノードとよばれます．

link tracing のイメージを図5に示しました．図では，s1, s2 の2つのソースノード，t1, t2 の2つのターゲットノードがまず抽出され，続いてソースからターゲットまでの最短経路にあるエッジとノードを抽出しています．

f:id:kur0cky:20201015112500p:plain — 図5：link tracing

link tracing が活躍するのは，ある現象の経路が簡単に手に入る場合です．たとえば，後ほど紹介する物理的なインターネットがその最たる例です．一つ一つのルーターのルーティング機能は単純で，通信経路は簡単なコマンドで知ることができますが，世界中に広がった全体の構造や性質はもはや個人には知ることができません．生き物などのエージェントの動きをトラッキングしている場合もそうでしょう．

ネットワークの大きさを推定する

前回記事では，母集団ネットワークの大きさ $N _ v$ は既知であるとしていました．たとえば日本人の交友関係に興味がある場合は，国勢調査などによりネットワークの大きさを予め知ることができます．そのため，ノードの包含確立 $\pi$ や二次の包含確率 $\pi _ {i, j}$ さえ分かれば，様々な量について Horvitz-Thompson 推定量を考えることができるのでした．しかし WWW など，ネットワークの大きさが未知であったり，それ自体が興味の対象となることも多いでしょう．

この節では，母集団ネットワークの大きさを推定する取り組みを2つ紹介します．これらは Kolaczyk(2009) でも引用されている取り組みですが，せっかくなのでもう少し詳しく紹介したいと思います．

link tracing を用いた物理インターネットの大きさの推定 (Viger et al., 2007)

Viger et al. (2007) は，特にルータなどによって構成される物理的なインターネットの大きさを推定することを目的とし，link tracing に基づいた推定量を提案しています． link tracing を利用するモチベーションに，Linux のtracerouteコマンドがあります． tracerouteでは，指定したホストに到達するまでの経路を取得することができます．コマンドを実行するIPをソースノードとし，指定したホストをターゲットノードとみなすと，これは link tracing そのものだとわかります．つまり，この論文では，traceroute を様々な端末から様々なホストに対して実行して得た標本ネットワークから，物理インターネットの大きさを推定することを目的としているのです．

論文で提案されている推定量は2つあり，Re Sampling (RS) estimator と Leave-One-Out (LOO) estimator です．どちらの推定量も，「ノードが重複して抽出される割合を用い，観測されたノードの数を膨らませる」というアイデアに基づいています．

Re-Sampling estimator

RS estimator では，標本ネットワーク $\mathcal G ^ {*}$ からさらに標本抽出を行い $\mathcal G ^ {**}$ を得ることが， $\mathcal G$ から $\mathcal G ^ {*}$ が観測されることと同様の性質を持っていると仮定します．より具体的には，ネットワーク $\mathcal G$ の大きさを $N _ {v}$ とすると， $N _ v$ と $N ^ * _ v$ ， $N ^ * _ v$ と $N ^ {* *} _ v$ の関係が等しいと仮定し，

$\displaystyle \begin{align} & \frac{E[N ^ {**} _ v ]}{N ^ {*} _ v} = \frac{N ^ {*} _ v}{N _ v}, \\\ \Leftrightarrow &N = N ^ {*} _ v \frac{N ^ {*} _ v}{E[N ^ {**} _ v ]}, \end{align}$

と整理できるため，次式で推定してやろうということです．

$\displaystyle \hat N _ {RS} = N ^ {*} \frac{N ^ {*}}{\bar N ^ {**}}.$

RS estimator で置かれている仮定はかなり厳しいものであると思います． $\mathcal G ^ *$ が $\mathcal G$ と同じ性質の構造になっているとは考えにくいでしょう．スケールフリーネットワークの部分グラフはスケールフリーにならない，という研究が Stumpf et al. (2005)などで行われていたりします．ましてや link tracing では，ハブとなるようなノードが抽出されやすく，得られる部分グラフは偏ったものと直感的に考えられます．

Leave-One-Out estimator

LOO estimator では，ターゲットノードが無作為抽出されていることを仮定します．そしてその名の通り，1つのターゲットノードを抜き出し，それが他のターゲットノードへの link tracing に含まれている割合を算出し，利用します，観測されたノードの数を膨らませます．

いま，ターゲットノード $j \in V$ を除いた link tracing で得られるノード集合を $\mathcal V ^ * _ {(-j)}$ ， $j$ が $\mathcal V ^ * _ {(-j)}$ に含まれないかどうかの指示関数を $\delta _ j = I(t _ j \notin \mathcal V _ {(-j)} ^ {*})$ ，他の $\mathcal V ^ * _ {(-j)}$ に含まれないターゲットノードの数を $X = \sum\delta _ j$ で表すことにします．ターゲットノード $j$ はランダムに選ばれているため，他のターゲットノードへのlink tracingで観測される確率は次式のようになります．ここで， $n _ S$ はソースノードの数， $n _ T$ はターゲットノードの数， $N _ {(-j)} ^ *$ は $\mathcal V _ {(-j)} ^ *$ の大きさです．

$\displaystyle P(\delta _ j = 1 | \mathcal V ^ * _ {(-j)}) = \frac{N _ v - N ^ * _ {(-j)}}{N _ v - n _ S - n _ T + 1}.$

さらに和をとることで， $X$ の期待値を求めることができます．

$\displaystyle E(X) = \sum _ {j = 1} ^ {n _ T} P(\delta _ j = 1 | \mathcal V ^ * _ {(-j)}) = n _ T\frac{N - E(N ^ * _ {(-)})}{N - n _ S - n _ T + 1}.$

これを $N$ について整理すると次式のようになります．

$\displaystyle N = \frac{n _ T E(N ^ * _ {(-)}) - (n _ S - n _ T + 1)E(X) }{n _ T - E(X)}.$

このままでも良さそうですが，分母 $n _ T - E(X)$ が0になると推定がおかしくなってしまうので， $n _ T - E(X)$ も直接推定し利用することを考え，論文では最終的に次の推定量が提案されています．

$\displaystyle \hat N _ {LOO} \approx (n _ S + n _ T) + \frac{N ^ * - (n _ S + n _ T)}{1 - X/(n _ T + 1)}.$

数値実験と，物理インターネットの大きさの推定

論文では，提案した推定量の良さを測るための数値実験が行われています．性質が既知であるネットワーク生成モデルからたくさんのネットワークを適用し，そこからの link tracing により擬似的な標本ネットワークを得ます．標本ネットワークからもとの大きさを推定し， $\hat N / N$ の平均で評価を行っています．生成モデルには，有名なBAモデルとERモデルが使われています（BAモデルやERモデルについては次回以降の記事で紹介することになると思います）．その結果，ソースノードやターゲットノードの数が極端に少なくない限り，LOO estimator の性質が良いと結論付けられています．

たしかに，RS estimator はその仮定が怪しい感じもありますが，LOO estimator の方はターゲットノードさえ無作為抽出できる状況であれば良さそうですね．

最後に，物理インターネットの大きさの推定を実際に行っています．推定量には数値実験で評価の良かった $\hat N _ {LOO}$ を使用し，結果として約7200万の推定値が算出されています．実際に真の値を確かめることはできないのですが，より直感的な方法であるpingコマンドとの比較が行われています． pingは，指定したIPに疎通できるかどうかを調べるコマンドです．さらに，IPの数は最大でも $2 ^ 32$ 通りです．そこで，論文では無作為に3726773個の IP に pingを送信し，61246個の有効な応答 (1.64%) を得ています．よって，ping を使った方法では，物理インターネットの大きさは $2 ^ 32 \times 0.0164 = 70583737$ で，約7000万となり， $\hat N _ {LOO}$ と近い値が得られます.

このように，ターゲットノードさえ無作為抽出することさえできれば，link tracing で得られた標本ネットワークからは悪くない精度でその大きさが推定できそうです．

snowball sampling に基づいた，ネットワークの大きさの推定 (Frank and Snijders, 1994)

次に，もう少し古めの論文を紹介します．

Frank and Snijders (1994) では，snowball sampling によって抽出された標本ネットワークから，母集団ネットワークの大きさを推定する推定量をたくさん（なんと7個も）提案しています．それらは大きく

Model-based な推定量
Design-based な推定量

に分けられます．ここでは，そのうちいくつかを取り上げて紹介することにします．

Model-based estimation

Model-basedな推定量では，以下のような状況を仮定します．つまり，母集団にERモデルを仮定するという意味で model based になっています．

snowball sampling の seed ノードが確率 $\alpha$ で無作為抽出されている
母集団ネットワークがランダムグラフであり，全てのノード間で確率 $\beta$ でエッジが張られている

母集団ネットワークのノード集合を $V = \{1, \dots, N _ v\}$ とし，ノード $i$ が抽出されたかどうかの指示関数を $x _ i$ とします．また，エッジ集合を $E = \{e _ {i, j} | i, j \in V\}$ とし，エッジ $e _ {i, j}$ が抽出されたかどうかの指示関数を $y _ {i, j}$ とします．さらに，seed 集合 $V _ {(0)} ^ *$ の数を $N _ {v(0)} ^ *$ ，seed 間のエッジ数を $N _ {e(0, 0)} ^ *$ とすると， $V _ {(0)} ^ *$ で条件付けした $N _ {e(0, 0)} ^ *$ は二項分布 $Binomial(N _ {v(0)} ^ *(N _ {v(0)} ^ *-1), \beta)$ に従います．また，条件付けない場合の確率は次のようになります．

$\displaystyle P(N _ {e(0,0)} ^ *) = \sum _ {N _ {v(0)} ^ * = 0} ^ {N _ v} {} _ {N _ v} {\rm C} _ {N _ {v(0)} ^ *} \alpha ^ {N _ {v(0)} ^ *} (1 - \alpha) ^ {N _ v - N^*_v(0)} \times {} _ {N _ {v(0)} ^ *(N _ {v(0)} ^ * - 1)} {\rm C} _ {N _ {e(0,0)} ^ *} \beta ^ {N _ {e(0,0)} ^ *} (1 - \beta) ^ {N _ {v(0)} ^ *(N _ {v(0)} ^ *-1) - N _ {e(0,0)} ^ *}.$

よって， $N _ {e(0,0)} ^ *$ の条件付き期待値や期待値は次のようになります．

$\displaystyle \begin{align} &E(N _ {e(0,0)} ^ * | N _ {v(0)} ^ \*) = N _ {v(0)} ^ *(N _ {v(0)} ^ * - 1)\beta, \\\ &E(N _ {e(0, 0)} ^ *) = N _ v (N _ v - 1) \alpha ^ 2 \beta. \end{align}$

次に，seed $V _ {(0)} ^ *$ からwave1で観測されるノード $V _ {(1)} ^ *$ に張られるエッジの数を $N ^ * _ {e(0,1)}$ とすると，先ほどと同様 $V _ {(0)} ^ *$ を条件付けた $N ^ * _ {e(0,1)}$ の確率は二項分布 $Binomial(N _ {v(0)} ^ *(N _ v - N _ {v(0)} ^ *), \beta)$ に従います．また，で条件付けない確率は次のようになります．

$\displaystyle \sum _ {N _ {v(0)} ^ * = 0} ^ {N _ v} {} _ {N _ v} {\rm C} _ {N _ {v(0)} ^ *} \alpha ^ {N _ {v(0)} ^ *} (1 - \alpha) ^ {N _ v - N _ {v(0)} ^ *} \times {} _ {N _ {v(0)} ^ *(N _ v - N _ {v(0)} ^ *)} {\rm }C _ {N ^ * _ {e(0,1)}} \beta ^ {N ^ * _ {e(0,1)}} (1 - \beta) ^ {N _ {v(0)} ^ *(N _ v - N _ {v(0)} ^ *) - N ^ * _ {e(0,1)}}$

よって， $N _ {e(0,1)} ^ *$ の条件付き期待値や期待値は次のようになります．

$\displaystyle \begin{align} &E(N _ {e(0,1)} ^ *| N _ {v(0)} ^ *) = N _ {v(0)} ^ * (N _ v - N _ {v(0)} ^ *) \beta, \\\ &E(N _ {e(0,1)} ^ *) = N _ v (N _ v - 1)\alpha(1-\alpha)\beta. \end{align}$

ここまでの条件付き期待値による連立方程式を解き， $N _ {v(0)} ^ *$ を条件付けたときの $N _ v, \beta$ のモーメント推定量は次式のようになります．

$\displaystyle \begin{align} \hat \beta &= r/n(n-1), \\\ \hat N _ v &= [ nr + (n-1)s ]. \end{align}$

同様に， $N _ {v(0)} ^ *$ を条件付けなかった場合は， $\alpha, \beta, N _ v$ のモーメント推定量は次式のようになります．

$\displaystyle \begin{align} \hat \alpha &= N _ {e(0, 0)} ^ */(N _ {e(0, 0)} ^ * + N _ {e(0, 1)} ^ *), \\\ \hat \beta &= N _ {e(0, 0)} ^ *(N _ {e(0, 0)} ^ * + N _ {e(0, 1)} ^ *) / N _ {v(0)} ^ * [ (N _ {v(0)} ^ *-1)N _ {e(0, 0)} ^ * + N _ {v(0)} ^ *N _ {e(0, 1)} ^ * ], \\\ \hat N _ {v} &= N _ {v(0)} ^ *(N _ {e(0, 0)} ^ * + N _ {e(0, 1)} ^ *) / N _ {e(0, 0)} ^ *. \end{align}$

Design-based estimatior

母集団ネットワークにランダムグラフが仮定できる場合，つまり論文でいう Model-based estimator は非常に簡単でした．では，ランダムグラフが仮定できない場合はどうでしょうか． Design-based estimator では，seed が無作為抽出されることのみを仮定します．

Model-based estimatorで用意した数式に加えて，wave 1 で新たに観測されたノードの数を $N ^ * _ {v(1)}$ とします．また，母集団ネットワーク全体の弧（エッジの片割れ）の大きさを $w$ とします． $N ^ * _ {v(0)}, N ^ * _ {e(0,0)}, N ^ * _ {e(0,1)}, N ^ * _ {v(1)}$ の期待値は $x _ j$ や $y _ {ij}$ を考えることで簡単に導出でき，次式のようになります．ここで， $B _ j$ は，ノード $j$ と接続するノード集合 $B _ j = \{i \in V: y _ {ij} = 1\}$ で， $b _ j$ はその大きさです．

$\displaystyle \begin{align} E(N ^ * _ {v(0)}) &= E\left(\sum _ {i} x _ i\right) = N _ v \alpha,\\\ E(N ^ * _ {e(0,0)}) &= E\left(\sum _ {i \neq j} x _ i x _ j y _ {ij}\right) = (w-N _ v)\alpha ^ 2, \\\ E (N ^ * _ {e(0,1)}) &= E\left(\sum _ {i\neq j} x _ i (1 - x _ j) y _ {ij}\right) = (w-N _ v) \alpha(1-\alpha), \\\ E (N ^ * _ {v(1)}) &= E\left(\sum _ {j}(\max _ {i\in B _ j} (x _ i - x _ j) )\right) \\\ & = E\left(\sum _ j(1-x _ j - \min _ {i\in B _ j}(1-x _ i))\right)\\\ & = N _ v(1-\alpha) - \sum _ j(1 - \alpha) ^ {b _ j}. \end{align}$

上式からわかるように，統計量 $N ^ * _ {v(0)}, N ^ * _ {e(0,0)}, N ^ * _ {e(0,1)}$ を使用すると，Model-basedのアプローチで導出された推定量が得られることになります．さらに，論文では少なくとも1つの他のseedと接続しているseedの数 $k = \sum _ j x _ j\max _ {i\neq j}y _ {ij} x _ i$ を導入し，別の推定量も導出しています． $k$ の期待値は次式のようになります．

$\displaystyle E(k) = N _ v \alpha - \alpha \sum _ j (1-\alpha) ^ {b _ j - 1}.$

$E(k) = \alpha E(k + m)$ であることを利用すると，次式の推定量を得ることができます．

$\hat N _ v = N ^ * _ {v(0)}(k + N ^ * _ {v(1)}) / k.$

数値実験

Frank and Snijders (1994) では，他にも Model-based の最尤推定量や，seedの抽出確率 $\alpha$ が既知であるときの Horvitz-Thompson 推定量なども提案されています．いくつかの推定量に対しては標準誤差の評価もされています．

これら多くの推定量を良さ測るためにa数値シミュレーションが行われています．次の3つのモデルからネットワークを生成し，それに対し snowball sampling で標本ネットワークを得て， $N _ v$ を推定します．大きさは $N _ v = 100, 1000$ ，seed の抽出確率は $N _ v = 100$ で $\alpha = 0.2$ ， $N _ v = 1000$ で $\alpha = 0.06$ ，と設定されています．

定数モデル：全てのノードに対して，入次数が5であるような有向グラフ（接続はランダムに選択）
ERモデル：ノード間にエッジが張られる確率があらゆるペアで $5 / (N _ v)$ であるような，ランダムグラフ
2-ブロックモデル：半分のノードの間では確率 $12 / (N _ v-2)$ ，もう半分のノードの間では確率 $6 / (N _ v-2)$ ，これらのグループ間では確率 $1/v$ でエッジが張られるようなブロックモデル

論文は1994年に出版されているため，スケールフリーネットワーク（次回以降に紹介予定です）を用いた実験は行われていません．さらに $N _ v$ が小さいこともあり，結果についてはここでは取りあげませんが，気になる方は元論文を参照ください．

おわりに

連載第2回の本記事では，より複雑なネットワーク標本抽出として，snowball sampling と link tracing を紹介しました．また，前回では既知としていた，ネットワークの大きさ $N _ v$ を推定する取り組みを2本ほど紹介しました．これらの論文では，link tracing のターゲットノードやsnowball sampling の seed を無作為抽出することを仮定していました．しかし，社会調査を目的としていない観察データなどを扱っていると，この仮定は非常に厳しいものでしょう．実際に利用する際は，交絡因子を調整したりすることが望ましいと考えられます．

やや文量が多くなってしまいましたが今回の記事はここまでです．標本抽出まわりは一旦これで終わりにし，次回からはネットワークの統計モデルに入っていきたいと思います．

参考文献

Frank, O., & Snijders, T. (1994). Estimating the size of hidden populations using snowball sampling. Journal of Official Statistics, 10, 53-67.
Kolaczyk, E. D. (2009). Statistical Analysis of Network Data. New York: Springer.
Stumpf, M. P., Wiuf, C., & May, R. M. (2005). Subnets of scale-free networks are not scale-free: Sampling properties of networks. Proceedings of the National Academy of Sciences, 102(12), 4221-4224.
Viger, F., Barrat, A., Dall’Asta, L., Zhang, C. H., & Kolaczyk, E. D. (2007). What is the real size of a sampled network? the case of the Internet. Physical Review E, 75(5), 056111.

▼これまでの記事

第1回ネットワークデータと標本調査(1)