1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
|
\begin{Verbatim}[commandchars=\\\{\}]
\PY{k+kn}{import} \PY{n+nn}{numpy} \PY{k+kn}{as} \PY{n+nn}{np}
\PY{k+kn}{from} \PY{n+nn}{numpy.linalg} \PY{k+kn}{import} \PY{n}{norm}
\PY{k+kn}{import} \PY{n+nn}{numpy.random} \PY{k+kn}{as} \PY{n+nn}{nr}
\PY{k+kn}{from} \PY{n+nn}{scipy.optimize} \PY{k+kn}{import} \PY{n}{minimize}
\PY{k+kn}{import} \PY{n+nn}{matplotlib.pyplot} \PY{k+kn}{as} \PY{n+nn}{plt}
\PY{k+kn}{import} \PY{n+nn}{seaborn}
\PY{k+kn}{from} \PY{n+nn}{random} \PY{k+kn}{import} \PY{n}{random}\PY{p}{,} \PY{n}{randint}
\PY{n}{seaborn}\PY{o}{.}\PY{n}{set\PYZus{}style}\PY{p}{(}\PY{l+s}{\PYZdq{}}\PY{l+s}{white}\PY{l+s}{\PYZdq{}}\PY{p}{)}
\PY{k}{def} \PY{n+nf}{likelihood}\PY{p}{(}\PY{n}{p}\PY{p}{,} \PY{n}{x}\PY{p}{,} \PY{n}{y}\PY{p}{)}\PY{p}{:}
\PY{n}{a} \PY{o}{=} \PY{n}{np}\PY{o}{.}\PY{n}{dot}\PY{p}{(}\PY{n}{x}\PY{p}{,} \PY{n}{p}\PY{p}{)}
\PY{k}{return} \PY{n}{np}\PY{o}{.}\PY{n}{log}\PY{p}{(}\PY{l+m+mf}{1.} \PY{o}{\PYZhy{}} \PY{n}{np}\PY{o}{.}\PY{n}{exp}\PY{p}{(}\PY{o}{\PYZhy{}}\PY{n}{a}\PY{p}{[}\PY{n}{y}\PY{p}{]}\PY{p}{)}\PY{p}{)}\PY{o}{.}\PY{n}{sum}\PY{p}{(}\PY{p}{)} \PY{o}{\PYZhy{}} \PY{n}{a}\PY{p}{[}\PY{o}{\PYZti{}}\PY{n}{y}\PY{p}{]}\PY{o}{.}\PY{n}{sum}\PY{p}{(}\PY{p}{)}
\PY{k}{def} \PY{n+nf}{likelihood\PYZus{}gradient}\PY{p}{(}\PY{n}{p}\PY{p}{,} \PY{n}{x}\PY{p}{,} \PY{n}{y}\PY{p}{)}\PY{p}{:}
\PY{n}{a} \PY{o}{=} \PY{n}{np}\PY{o}{.}\PY{n}{dot}\PY{p}{(}\PY{n}{x}\PY{p}{,} \PY{n}{p}\PY{p}{)}
\PY{n}{l} \PY{o}{=} \PY{n}{np}\PY{o}{.}\PY{n}{log}\PY{p}{(}\PY{l+m+mf}{1.} \PY{o}{\PYZhy{}} \PY{n}{np}\PY{o}{.}\PY{n}{exp}\PY{p}{(}\PY{o}{\PYZhy{}}\PY{n}{a}\PY{p}{[}\PY{n}{y}\PY{p}{]}\PY{p}{)}\PY{p}{)}\PY{o}{.}\PY{n}{sum}\PY{p}{(}\PY{p}{)} \PY{o}{\PYZhy{}} \PY{n}{a}\PY{p}{[}\PY{o}{\PYZti{}}\PY{n}{y}\PY{p}{]}\PY{o}{.}\PY{n}{sum}\PY{p}{(}\PY{p}{)}
\PY{n}{g1} \PY{o}{=} \PY{l+m+mf}{1.} \PY{o}{/} \PY{p}{(}\PY{n}{np}\PY{o}{.}\PY{n}{exp}\PY{p}{(}\PY{n}{a}\PY{p}{[}\PY{n}{y}\PY{p}{]}\PY{p}{)} \PY{o}{\PYZhy{}} \PY{l+m+mf}{1.}\PY{p}{)}
\PY{n}{g} \PY{o}{=} \PY{p}{(}\PY{n}{x}\PY{p}{[}\PY{n}{y}\PY{p}{]} \PY{o}{*} \PY{n}{g1}\PY{p}{[}\PY{p}{:}\PY{p}{,} \PY{n}{np}\PY{o}{.}\PY{n}{newaxis}\PY{p}{]}\PY{p}{)}\PY{o}{.}\PY{n}{sum}\PY{p}{(}\PY{l+m+mi}{0}\PY{p}{)} \PY{o}{\PYZhy{}} \PY{n}{x}\PY{p}{[}\PY{o}{\PYZti{}}\PY{n}{y}\PY{p}{]}\PY{o}{.}\PY{n}{sum}\PY{p}{(}\PY{l+m+mi}{0}\PY{p}{)}
\PY{k}{return} \PY{n}{l}\PY{p}{,} \PY{n}{g}
\PY{k}{def} \PY{n+nf}{test\PYZus{}gradient}\PY{p}{(}\PY{n}{x}\PY{p}{,} \PY{n}{y}\PY{p}{)}\PY{p}{:}
\PY{n}{eps} \PY{o}{=} \PY{l+m+mf}{1e\PYZhy{}10}
\PY{k}{for} \PY{n}{i} \PY{o+ow}{in} \PY{n+nb}{xrange}\PY{p}{(}\PY{n}{x}\PY{o}{.}\PY{n}{shape}\PY{p}{[}\PY{l+m+mi}{1}\PY{p}{]}\PY{p}{)}\PY{p}{:}
\PY{n}{p} \PY{o}{=} \PY{l+m+mf}{0.5} \PY{o}{*} \PY{n}{np}\PY{o}{.}\PY{n}{ones}\PY{p}{(}\PY{n}{x}\PY{o}{.}\PY{n}{shape}\PY{p}{[}\PY{l+m+mi}{1}\PY{p}{]}\PY{p}{)}
\PY{n}{a} \PY{o}{=} \PY{n}{np}\PY{o}{.}\PY{n}{dot}\PY{p}{(}\PY{n}{x}\PY{p}{,} \PY{n}{p}\PY{p}{)}
\PY{n}{g1} \PY{o}{=} \PY{l+m+mf}{1.} \PY{o}{/} \PY{p}{(}\PY{n}{np}\PY{o}{.}\PY{n}{exp}\PY{p}{(}\PY{n}{a}\PY{p}{[}\PY{n}{y}\PY{p}{]}\PY{p}{)} \PY{o}{\PYZhy{}} \PY{l+m+mf}{1.}\PY{p}{)}
\PY{n}{g} \PY{o}{=} \PY{p}{(}\PY{n}{x}\PY{p}{[}\PY{n}{y}\PY{p}{]} \PY{o}{*} \PY{n}{g1}\PY{p}{[}\PY{p}{:}\PY{p}{,} \PY{n}{np}\PY{o}{.}\PY{n}{newaxis}\PY{p}{]}\PY{p}{)}\PY{o}{.}\PY{n}{sum}\PY{p}{(}\PY{l+m+mi}{0}\PY{p}{)} \PY{o}{\PYZhy{}} \PY{n}{x}\PY{p}{[}\PY{o}{\PYZti{}}\PY{n}{y}\PY{p}{]}\PY{o}{.}\PY{n}{sum}\PY{p}{(}\PY{l+m+mi}{0}\PY{p}{)}
\PY{n}{p}\PY{p}{[}\PY{n}{i}\PY{p}{]} \PY{o}{+}\PY{o}{=} \PY{n}{eps}
\PY{n}{f1} \PY{o}{=} \PY{n}{likelihood}\PY{p}{(}\PY{n}{p}\PY{p}{,} \PY{n}{x}\PY{p}{,} \PY{n}{y}\PY{p}{)}
\PY{n}{p}\PY{p}{[}\PY{n}{i}\PY{p}{]} \PY{o}{\PYZhy{}}\PY{o}{=} \PY{l+m+mi}{2} \PY{o}{*} \PY{n}{eps}
\PY{n}{f2} \PY{o}{=} \PY{n}{likelihood}\PY{p}{(}\PY{n}{p}\PY{p}{,} \PY{n}{x}\PY{p}{,} \PY{n}{y}\PY{p}{)}
\PY{k}{print} \PY{n}{g}\PY{p}{[}\PY{n}{i}\PY{p}{]}\PY{p}{,} \PY{p}{(}\PY{n}{f1} \PY{o}{\PYZhy{}} \PY{n}{f2}\PY{p}{)} \PY{o}{/} \PY{p}{(}\PY{l+m+mi}{2} \PY{o}{*} \PY{n}{eps}\PY{p}{)}
\PY{k}{def} \PY{n+nf}{infer}\PY{p}{(}\PY{n}{x}\PY{p}{,} \PY{n}{y}\PY{p}{)}\PY{p}{:}
\PY{k}{def} \PY{n+nf}{f}\PY{p}{(}\PY{n}{p}\PY{p}{)}\PY{p}{:}
\PY{n}{l}\PY{p}{,} \PY{n}{g} \PY{o}{=} \PY{n}{likelihood\PYZus{}gradient}\PY{p}{(}\PY{n}{p}\PY{p}{,} \PY{n}{x}\PY{p}{,} \PY{n}{y}\PY{p}{)}
\PY{k}{return} \PY{o}{\PYZhy{}}\PY{n}{l}\PY{p}{,} \PY{o}{\PYZhy{}}\PY{n}{g}
\PY{n}{x0} \PY{o}{=} \PY{n}{np}\PY{o}{.}\PY{n}{ones}\PY{p}{(}\PY{n}{x}\PY{o}{.}\PY{n}{shape}\PY{p}{[}\PY{l+m+mi}{1}\PY{p}{]}\PY{p}{)}
\PY{n}{bounds} \PY{o}{=} \PY{p}{[}\PY{p}{(}\PY{l+m+mf}{1e\PYZhy{}10}\PY{p}{,} \PY{n+nb+bp}{None}\PY{p}{)}\PY{p}{]} \PY{o}{*} \PY{n}{x}\PY{o}{.}\PY{n}{shape}\PY{p}{[}\PY{l+m+mi}{1}\PY{p}{]}
\PY{k}{return} \PY{n}{minimize}\PY{p}{(}\PY{n}{f}\PY{p}{,} \PY{n}{x0}\PY{p}{,} \PY{n}{jac}\PY{o}{=}\PY{n+nb+bp}{True}\PY{p}{,} \PY{n}{bounds}\PY{o}{=}\PY{n}{bounds}\PY{p}{,} \PY{n}{method}\PY{o}{=}\PY{l+s}{\PYZdq{}}\PY{l+s}{L\PYZhy{}BFGS\PYZhy{}B}\PY{l+s}{\PYZdq{}}\PY{p}{)}\PY{o}{.}\PY{n}{x}
\PY{k}{def} \PY{n+nf}{bootstrap}\PY{p}{(}\PY{n}{x}\PY{p}{,} \PY{n}{y}\PY{p}{,} \PY{n}{n\PYZus{}iter}\PY{o}{=}\PY{l+m+mi}{100}\PY{p}{)}\PY{p}{:}
\PY{n}{rval} \PY{o}{=} \PY{n}{np}\PY{o}{.}\PY{n}{zeros}\PY{p}{(}\PY{p}{(}\PY{n}{n\PYZus{}iter}\PY{p}{,} \PY{n}{x}\PY{o}{.}\PY{n}{shape}\PY{p}{[}\PY{l+m+mi}{1}\PY{p}{]}\PY{p}{)}\PY{p}{)}
\PY{k}{for} \PY{n}{i} \PY{o+ow}{in} \PY{n+nb}{xrange}\PY{p}{(}\PY{n}{n\PYZus{}iter}\PY{p}{)}\PY{p}{:}
\PY{n}{indices} \PY{o}{=} \PY{n}{np}\PY{o}{.}\PY{n}{random}\PY{o}{.}\PY{n}{choice}\PY{p}{(}\PY{n+nb}{len}\PY{p}{(}\PY{n}{y}\PY{p}{)}\PY{p}{,} \PY{n}{replace}\PY{o}{=}\PY{n+nb+bp}{False}\PY{p}{,} \PY{n}{size}\PY{o}{=}\PY{n+nb}{int}\PY{p}{(}\PY{n+nb}{len}\PY{p}{(}\PY{n}{y}\PY{p}{)}\PY{o}{*}\PY{o}{.}\PY{l+m+mi}{9}\PY{p}{)}\PY{p}{)}
\PY{n}{rval}\PY{p}{[}\PY{n}{i}\PY{p}{]} \PY{o}{=} \PY{n}{infer}\PY{p}{(}\PY{n}{x}\PY{p}{[}\PY{n}{indices}\PY{p}{]}\PY{p}{,} \PY{n}{y}\PY{p}{[}\PY{n}{indices}\PY{p}{]}\PY{p}{)}
\PY{k}{return} \PY{n}{rval}
\PY{k}{def} \PY{n+nf}{confidence\PYZus{}interval}\PY{p}{(}\PY{n}{counts}\PY{p}{,} \PY{n}{bins}\PY{p}{)}\PY{p}{:}
\PY{n}{k} \PY{o}{=} \PY{l+m+mi}{0}
\PY{k}{while} \PY{n}{np}\PY{o}{.}\PY{n}{sum}\PY{p}{(}\PY{n}{counts}\PY{p}{[}\PY{n+nb}{len}\PY{p}{(}\PY{n}{counts}\PY{p}{)}\PY{o}{/}\PY{l+m+mi}{2}\PY{o}{\PYZhy{}}\PY{n}{k}\PY{p}{:}\PY{n+nb}{len}\PY{p}{(}\PY{n}{counts}\PY{p}{)}\PY{o}{/}\PY{l+m+mi}{2}\PY{o}{+}\PY{n}{k}\PY{p}{]}\PY{p}{)} \PY{o}{\PYZlt{}}\PY{o}{=} \PY{o}{.}\PY{l+m+mi}{95}\PY{o}{*}\PY{n}{np}\PY{o}{.}\PY{n}{sum}\PY{p}{(}\PY{n}{counts}\PY{p}{)}\PY{p}{:}
\PY{n}{k} \PY{o}{+}\PY{o}{=} \PY{l+m+mi}{1}
\PY{k}{return} \PY{n}{bins}\PY{p}{[}\PY{n+nb}{len}\PY{p}{(}\PY{n}{bins}\PY{p}{)}\PY{o}{/}\PY{l+m+mi}{2}\PY{o}{\PYZhy{}}\PY{n}{k}\PY{p}{]}\PY{p}{,} \PY{n}{bins}\PY{p}{[}\PY{n+nb}{len}\PY{p}{(}\PY{n}{bins}\PY{p}{)}\PY{o}{/}\PY{l+m+mi}{2}\PY{o}{+}\PY{n}{k}\PY{p}{]}
\PY{k}{def} \PY{n+nf}{build\PYZus{}matrix}\PY{p}{(}\PY{n}{cascades}\PY{p}{,} \PY{n}{node}\PY{p}{)}\PY{p}{:}
\PY{k}{def} \PY{n+nf}{aux}\PY{p}{(}\PY{n}{cascade}\PY{p}{,} \PY{n}{node}\PY{p}{)}\PY{p}{:}
\PY{n}{xlist}\PY{p}{,} \PY{n}{slist} \PY{o}{=} \PY{n+nb}{zip}\PY{p}{(}\PY{o}{*}\PY{n}{cascade}\PY{p}{)}
\PY{n}{indices} \PY{o}{=} \PY{p}{[}\PY{n}{i} \PY{k}{for} \PY{n}{i}\PY{p}{,} \PY{n}{s} \PY{o+ow}{in} \PY{n+nb}{enumerate}\PY{p}{(}\PY{n}{slist}\PY{p}{)} \PY{k}{if} \PY{n}{s}\PY{p}{[}\PY{n}{node}\PY{p}{]} \PY{o+ow}{and} \PY{n}{i} \PY{o}{\PYZgt{}}\PY{o}{=} \PY{l+m+mi}{1}\PY{p}{]}
\PY{k}{if} \PY{n}{indices}\PY{p}{:}
\PY{n}{x} \PY{o}{=} \PY{n}{np}\PY{o}{.}\PY{n}{vstack}\PY{p}{(}\PY{n}{xlist}\PY{p}{[}\PY{n}{i}\PY{o}{\PYZhy{}}\PY{l+m+mi}{1}\PY{p}{]} \PY{k}{for} \PY{n}{i} \PY{o+ow}{in} \PY{n}{indices}\PY{p}{)}
\PY{n}{y} \PY{o}{=} \PY{n}{np}\PY{o}{.}\PY{n}{array}\PY{p}{(}\PY{p}{[}\PY{n}{xlist}\PY{p}{[}\PY{n}{i}\PY{p}{]}\PY{p}{[}\PY{n}{node}\PY{p}{]} \PY{k}{for} \PY{n}{i} \PY{o+ow}{in} \PY{n}{indices}\PY{p}{]}\PY{p}{)}
\PY{k}{return} \PY{n}{x}\PY{p}{,} \PY{n}{y}
\PY{k}{else}\PY{p}{:}
\PY{k}{return} \PY{n+nb+bp}{None}
\PY{n}{pairs} \PY{o}{=} \PY{p}{(}\PY{n}{aux}\PY{p}{(}\PY{n}{cascade}\PY{p}{,} \PY{n}{node}\PY{p}{)} \PY{k}{for} \PY{n}{cascade} \PY{o+ow}{in} \PY{n}{cascades}\PY{p}{)}
\PY{n}{xs}\PY{p}{,} \PY{n}{ys} \PY{o}{=} \PY{n+nb}{zip}\PY{p}{(}\PY{o}{*}\PY{p}{(}\PY{n}{pair} \PY{k}{for} \PY{n}{pair} \PY{o+ow}{in} \PY{n}{pairs} \PY{k}{if} \PY{n}{pair}\PY{p}{)}\PY{p}{)}
\PY{n}{x} \PY{o}{=} \PY{n}{np}\PY{o}{.}\PY{n}{vstack}\PY{p}{(}\PY{n}{xs}\PY{p}{)}
\PY{n}{y} \PY{o}{=} \PY{n}{np}\PY{o}{.}\PY{n}{concatenate}\PY{p}{(}\PY{n}{ys}\PY{p}{)}
\PY{k}{return} \PY{n}{x}\PY{p}{,} \PY{n}{y}
\PY{k}{def} \PY{n+nf}{build\PYZus{}cascade\PYZus{}list}\PY{p}{(}\PY{n}{cascades}\PY{p}{)}\PY{p}{:}
\PY{n}{x}\PY{p}{,} \PY{n}{s} \PY{o}{=} \PY{p}{[}\PY{p}{]}\PY{p}{,} \PY{p}{[}\PY{p}{]}
\PY{k}{for} \PY{n}{cascade} \PY{o+ow}{in} \PY{n}{cascades}\PY{p}{:}
\PY{n}{xlist}\PY{p}{,} \PY{n}{slist} \PY{o}{=} \PY{n+nb}{zip}\PY{p}{(}\PY{o}{*}\PY{n}{cascade}\PY{p}{)}
\PY{n}{x}\PY{o}{.}\PY{n}{append}\PY{p}{(}\PY{n}{xlist}\PY{p}{)}
\PY{n}{s}\PY{o}{.}\PY{n}{append}\PY{p}{(}\PY{n}{slist}\PY{p}{)}
\PY{k}{return} \PY{n}{x}\PY{p}{,} \PY{n}{s}
\PY{k}{def} \PY{n+nf}{simulate\PYZus{}cascade}\PY{p}{(}\PY{n}{x}\PY{p}{,} \PY{n}{graph}\PY{p}{)}\PY{p}{:}
\PY{l+s+sd}{\PYZdq{}\PYZdq{}\PYZdq{}}
\PY{l+s+sd}{ Simulate an IC cascade given a graph and initial state.}
\PY{l+s+sd}{ For each time step we yield:}
\PY{l+s+sd}{ \PYZhy{} susc: the nodes susceptible at the beginning of this time step}
\PY{l+s+sd}{ \PYZhy{} x: the subset of susc who became infected}
\PY{l+s+sd}{ \PYZdq{}\PYZdq{}\PYZdq{}}
\PY{n}{susc} \PY{o}{=} \PY{n}{np}\PY{o}{.}\PY{n}{ones}\PY{p}{(}\PY{n}{graph}\PY{o}{.}\PY{n}{shape}\PY{p}{[}\PY{l+m+mi}{0}\PY{p}{]}\PY{p}{,} \PY{n}{dtype}\PY{o}{=}\PY{n+nb}{bool}\PY{p}{)} \PY{c}{\PYZsh{} t=0, everyone is susceptible}
\PY{k}{yield} \PY{n}{x}\PY{p}{,} \PY{n}{susc}
\PY{k}{while} \PY{n}{np}\PY{o}{.}\PY{n}{any}\PY{p}{(}\PY{n}{x}\PY{p}{)}\PY{p}{:}
\PY{n}{susc} \PY{o}{=} \PY{n}{susc} \PY{o}{\PYZca{}} \PY{n}{x} \PY{c}{\PYZsh{} nodes infected at previous step are now inactive}
\PY{k}{if} \PY{o+ow}{not} \PY{n}{np}\PY{o}{.}\PY{n}{any}\PY{p}{(}\PY{n}{susc}\PY{p}{)}\PY{p}{:}
\PY{k}{break}
\PY{n}{x} \PY{o}{=} \PY{l+m+mi}{1} \PY{o}{\PYZhy{}} \PY{n}{np}\PY{o}{.}\PY{n}{exp}\PY{p}{(}\PY{o}{\PYZhy{}}\PY{n}{np}\PY{o}{.}\PY{n}{dot}\PY{p}{(}\PY{n}{graph}\PY{o}{.}\PY{n}{T}\PY{p}{,} \PY{n}{x}\PY{p}{)}\PY{p}{)}
\PY{n}{y} \PY{o}{=} \PY{n}{nr}\PY{o}{.}\PY{n}{random}\PY{p}{(}\PY{n}{x}\PY{o}{.}\PY{n}{shape}\PY{p}{[}\PY{l+m+mi}{0}\PY{p}{]}\PY{p}{)}
\PY{n}{x} \PY{o}{=} \PY{p}{(}\PY{n}{x} \PY{o}{\PYZgt{}}\PY{o}{=} \PY{n}{y}\PY{p}{)} \PY{o}{\PYZam{}} \PY{n}{susc}
\PY{k}{yield} \PY{n}{x}\PY{p}{,} \PY{n}{susc}
\PY{k}{def} \PY{n+nf}{uniform\PYZus{}source}\PY{p}{(}\PY{n}{graph}\PY{p}{)}\PY{p}{:}
\PY{n}{x0} \PY{o}{=} \PY{n}{np}\PY{o}{.}\PY{n}{zeros}\PY{p}{(}\PY{n}{graph}\PY{o}{.}\PY{n}{shape}\PY{p}{[}\PY{l+m+mi}{0}\PY{p}{]}\PY{p}{,} \PY{n}{dtype}\PY{o}{=}\PY{n+nb}{bool}\PY{p}{)}
\PY{n}{x0}\PY{p}{[}\PY{n}{nr}\PY{o}{.}\PY{n}{randint}\PY{p}{(}\PY{l+m+mi}{0}\PY{p}{,} \PY{n}{graph}\PY{o}{.}\PY{n}{shape}\PY{p}{[}\PY{l+m+mi}{0}\PY{p}{]}\PY{p}{)}\PY{p}{]} \PY{o}{=} \PY{n+nb+bp}{True}
\PY{k}{return} \PY{n}{x0}
\PY{k}{def} \PY{n+nf}{simulate\PYZus{}cascades}\PY{p}{(}\PY{n}{n}\PY{p}{,} \PY{n}{graph}\PY{p}{,} \PY{n}{source}\PY{o}{=}\PY{n}{uniform\PYZus{}source}\PY{p}{)}\PY{p}{:}
\PY{k}{for} \PY{n}{\PYZus{}} \PY{o+ow}{in} \PY{n+nb}{xrange}\PY{p}{(}\PY{n}{n}\PY{p}{)}\PY{p}{:}
\PY{n}{x0} \PY{o}{=} \PY{n}{source}\PY{p}{(}\PY{n}{graph}\PY{p}{)}
\PY{k}{yield} \PY{n}{simulate\PYZus{}cascade}\PY{p}{(}\PY{n}{x0}\PY{p}{,} \PY{n}{graph}\PY{p}{)}
\PY{k}{if} \PY{n}{\PYZus{}\PYZus{}name\PYZus{}\PYZus{}} \PY{o}{==} \PY{l+s}{\PYZdq{}}\PY{l+s}{\PYZus{}\PYZus{}main\PYZus{}\PYZus{}}\PY{l+s}{\PYZdq{}}\PY{p}{:}
\PY{c}{\PYZsh{} g = np.array([[0, 1, 1, 0], [1, 0, 0, 1], [1, 0, 0, 1], [0, 1, 1, 0]])}
\PY{n}{g} \PY{o}{=} \PY{n}{np}\PY{o}{.}\PY{n}{array}\PY{p}{(}\PY{p}{[}\PY{p}{[}\PY{l+m+mi}{0}\PY{p}{,} \PY{l+m+mi}{0}\PY{p}{,} \PY{l+m+mi}{1}\PY{p}{]}\PY{p}{,} \PY{p}{[}\PY{l+m+mi}{0}\PY{p}{,} \PY{l+m+mi}{0}\PY{p}{,} \PY{l+m+mf}{0.5}\PY{p}{]}\PY{p}{,} \PY{p}{[}\PY{l+m+mi}{0}\PY{p}{,} \PY{l+m+mi}{0}\PY{p}{,} \PY{l+m+mi}{0}\PY{p}{]}\PY{p}{]}\PY{p}{)}
\PY{n}{p} \PY{o}{=} \PY{l+m+mf}{0.5}
\PY{n}{g} \PY{o}{=} \PY{n}{np}\PY{o}{.}\PY{n}{log}\PY{p}{(}\PY{l+m+mf}{1.} \PY{o}{/} \PY{p}{(}\PY{l+m+mi}{1} \PY{o}{\PYZhy{}} \PY{n}{p} \PY{o}{*} \PY{n}{g}\PY{p}{)}\PY{p}{)}
\PY{c}{\PYZsh{} error = []}
\PY{k}{def} \PY{n+nf}{source}\PY{p}{(}\PY{n}{graph}\PY{p}{,} \PY{n}{t}\PY{p}{)}\PY{p}{:}
\PY{n}{x0} \PY{o}{=} \PY{n}{np}\PY{o}{.}\PY{n}{zeros}\PY{p}{(}\PY{n}{graph}\PY{o}{.}\PY{n}{shape}\PY{p}{[}\PY{l+m+mi}{0}\PY{p}{]}\PY{p}{,} \PY{n}{dtype}\PY{o}{=}\PY{n+nb}{bool}\PY{p}{)}
\PY{n}{a} \PY{o}{=} \PY{n}{randint}\PY{p}{(}\PY{l+m+mi}{0}\PY{p}{,} \PY{l+m+mi}{1}\PY{p}{)}
\PY{n}{x0}\PY{p}{[}\PY{n}{a}\PY{p}{]} \PY{o}{=} \PY{n+nb+bp}{True}
\PY{k}{if} \PY{n}{random}\PY{p}{(}\PY{p}{)} \PY{o}{\PYZgt{}} \PY{n}{t}\PY{p}{:}
\PY{n}{x0}\PY{p}{[}\PY{l+m+mi}{1}\PY{o}{\PYZhy{}}\PY{n}{a}\PY{p}{]} \PY{o}{=} \PY{n+nb+bp}{True}
\PY{k}{return} \PY{n}{x0}
\PY{n}{thresh} \PY{o}{=} \PY{n}{np}\PY{o}{.}\PY{n}{arange}\PY{p}{(}\PY{l+m+mf}{0.}\PY{p}{,} \PY{l+m+mf}{1.1}\PY{p}{,} \PY{n}{step}\PY{o}{=}\PY{l+m+mf}{0.2}\PY{p}{)}
\PY{n}{sizes} \PY{o}{=} \PY{n}{np}\PY{o}{.}\PY{n}{arange}\PY{p}{(}\PY{l+m+mi}{10}\PY{p}{,} \PY{l+m+mi}{100}\PY{p}{,} \PY{n}{step}\PY{o}{=}\PY{l+m+mi}{10}\PY{p}{)}
\PY{n}{nsimul} \PY{o}{=} \PY{l+m+mi}{10}
\PY{n}{r} \PY{o}{=} \PY{n}{np}\PY{o}{.}\PY{n}{zeros}\PY{p}{(}\PY{n+nb}{len}\PY{p}{(}\PY{n}{sizes}\PY{p}{)}\PY{p}{,} \PY{n+nb}{len}\PY{p}{(}\PY{n}{thresh}\PY{p}{)}\PY{p}{)}
\PY{k}{for} \PY{n}{t} \PY{o+ow}{in} \PY{n}{thresh}\PY{p}{:}
\PY{k}{for} \PY{n}{i} \PY{o+ow}{in} \PY{n}{nsimul}\PY{p}{:}
\PY{n}{cascades} \PY{o}{=} \PY{n}{simulate\PYZus{}cascades}\PY{p}{(}\PY{n}{np}\PY{o}{.}\PY{n}{max}\PY{p}{(}\PY{n}{sizes}\PY{p}{)}\PY{p}{,} \PY{n}{g}\PY{p}{,}
\PY{n}{source}\PY{o}{=}\PY{k}{lambda} \PY{n}{graph}\PY{p}{:} \PY{n}{source}\PY{p}{(}\PY{n}{graph}\PY{p}{,} \PY{n}{t}\PY{p}{)}\PY{p}{)}
\PY{n}{e} \PY{o}{=} \PY{n}{np}\PY{o}{.}\PY{n}{zeros}\PY{p}{(}\PY{n}{g}\PY{o}{.}\PY{n}{shape}\PY{p}{[}\PY{l+m+mi}{0}\PY{p}{]}\PY{p}{)}
\PY{k}{for} \PY{n}{j}\PY{p}{,} \PY{n}{s} \PY{o+ow}{in} \PY{n+nb}{enumerate}\PY{p}{(}\PY{n}{sizes}\PY{p}{)}\PY{p}{:}
\PY{n}{x}\PY{p}{,} \PY{n}{y} \PY{o}{=} \PY{n}{build\PYZus{}matrix}\PY{p}{(}\PY{n}{cascades}\PY{p}{,} \PY{l+m+mi}{2}\PY{p}{)}
\PY{n}{e} \PY{o}{+}\PY{o}{=} \PY{n}{infer}\PY{p}{(}\PY{n}{x}\PY{p}{[}\PY{p}{:}\PY{n}{s}\PY{p}{]}\PY{p}{,} \PY{n}{y}\PY{p}{[}\PY{p}{:}\PY{n}{s}\PY{p}{]}\PY{p}{)}
\PY{k}{for} \PY{n}{i}\PY{p}{,} \PY{n}{t} \PY{o+ow}{in} \PY{n+nb}{enumerate}\PY{p}{(}\PY{n}{thresh}\PY{p}{)}\PY{p}{:}
\PY{n}{plt}\PY{o}{.}\PY{n}{plot}\PY{p}{(}\PY{n}{sizes}\PY{p}{,} \PY{n}{m}\PY{p}{[}\PY{p}{:}\PY{p}{,} \PY{n}{i}\PY{p}{]}\PY{p}{,} \PY{n}{label}\PY{o}{=}\PY{n+nb}{str}\PY{p}{(}\PY{n}{t}\PY{p}{)}\PY{p}{)}
\PY{n}{plt}\PY{o}{.}\PY{n}{legend}\PY{p}{(}\PY{p}{)}
\PY{n}{plt}\PY{o}{.}\PY{n}{show}\PY{p}{(}\PY{p}{)}
\PY{c}{\PYZsh{} conf = bootstrap(x, y, n\PYZus{}iter=100)}
\PY{c}{\PYZsh{} estimand = np.linalg.norm(np.delete(conf \PYZhy{} g[0], 0, axis=1), axis=1)}
\PY{c}{\PYZsh{} error.append(confidence\PYZus{}interval(*np.histogram(estimand, bins=50)))}
\PY{c}{\PYZsh{} plt.semilogx(sizes, error)}
\PY{c}{\PYZsh{} plt.show()}
\end{Verbatim}
|
