aboutsummaryrefslogtreecommitdiffstats
path: root/finale/mle.tex
diff options
context:
space:
mode:
Diffstat (limited to 'finale/mle.tex')
-rw-r--r--finale/mle.tex163
1 files changed, 163 insertions, 0 deletions
diff --git a/finale/mle.tex b/finale/mle.tex
new file mode 100644
index 0000000..70caa3e
--- /dev/null
+++ b/finale/mle.tex
@@ -0,0 +1,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}