\documentclass[10pt]{article}
\usepackage{fullpage, amsmath, amsthm, amssymb}
\usepackage{graphicx, algorithm, algpseudocode}
\usepackage[utf8x]{inputenc}
\newtheorem{theorem}{Theorem}

\title{Sparse Recovery for Graph Reconstruction}
\author{Eric Balkanski \and Jean Pouget-Abadie}

\date{}

\begin{document}

\maketitle

\section{Introduction}

Given a set of observed cascades, the graph reconstruction problem consists of finding the underlying graph on which these cascades spread. We explore different models for cascade creation, namely the voter model and the independent cascade model. Our objective is therefore to approximately compute the `inverse' of the cascade generation function from as few observations as possible.

\section{Related Work}

There have been several works tackling the graph reconstruction problem in variants of the independent cascade. We briefly summarize their results and approaches below.

\begin{itemize}
\item Manuel Gomez-Rodriguez, Jure Leskovec, and Andreas Klause were amongst the first to tackle the problem of graph reconstruction from cascade observation \cite{gomez}. Their first method \textsc{NetInf}, and a series of other similar methods published subsequently, solves a heuristic which approximates a NP-hard maximum-likelihood formulation. Unfortunately \textsc{NetInf} has no known theoretical guarantees. One difficulty with comparing ourselves to their heuristic is that they assume a continuous time independent cascade model, which is close but un-adaptable to our framework. Namely, as we will see in Section~\ref{sec:icc}, nodes in our model cannot influence nodes beyond the time step at which they are active.
\item Their paper was later followed by a publication from Praneeth Netrapalli and Sujay Sanghavi \cite{netrapalli}, which provided two algorithms with several theoretical guarantees. 
  \begin{itemize}
  \item The first algorithm they discuss is a convex and parallelizable optimization problem. However, the theoretical guarantee is close to perfect reconstruction of `strong' edges after ${\cal O}(d^2 \log n)$ with high probability. There are several problems to this result however.
  \begin{itemize}
    \item Researchers believe the true threshold of cascades necessary is closer to ${\cal O}(d \log n)$
    \item The constant in the big ${\cal O}$ is highly polynomial in the different constants, of the order of billions even for very generous initializations of these parameters
    \item The authors make a restrictive assumption on the probabilities of each edge, which they call `correlation decay': $\forall i, \sum_{j \in {\cal N}(i)} p_{j,i} < 1 -\alpha$, where $0 < \alpha < 1$
  \end{itemize}
  \item The second algorithm, to which we will compare ourselves in Section~\ref{sec:exp}, is \textsc{Greedy}. It achieves perfect reconstruction after ${\cal O}(d \log n)$ observed cascades with high probability. However, this guarantee has been proven only in the case of tree-graphs. The authors mention that it does fairly well in pratice on other types of graphs.
  \end{itemize}
\item Finally, a recent paper by Bruno Abrahao, Flavio Chierichetti, Robert Kleinberg, Alessandro Panconesi gives several algorithms and a different approach to showing their theoretical guarantees. Their strongest result is for bounded-degree graphs, where they show that they can obtain perfect reconstruction after observing ${\cal O}(\Delta^9 \log^2 \Delta \log n)$ algorithm. This is weaker than the first algorithm suggested by \cite{netrapalli}, but does not make the `correlation decay' assumption., where $\Delta$ is the max degree in graph ${\cal G}$.
\end{itemize}

What we aim for with the following framework is to achieve theoretical guarantees promising close-to-perfect reconstruction with observing ${\cal O}(\Delta \log n)$ cascades with high probability without relying on unreasonable assumptions like `correlation decay'. We also want our method to be easily implementable in practice and compare favorably to the above algorithms.

\section{The Voter Model}

\subsection{Description}
In the voter model, there are two types of nodes, \textsc{red} and \textsc{blue}. We fix a horizon $T > 0$. At $t=0$, we decide on a set of \textsc{blue} nodes, called the seed set. Every other node is colored $\textsc{red}$. At each time step, each node $u$ chooses one of its neighbors uniformly, with probability $1/deg(u)$, and adopts the color of that neighbor. A cascade in the voter model can be described as the evolution of which nodes are \textsc{blue} for each time step $0\leq t\leq T$.

Several variants of the model can be formulated, namely whether or not nodes have self-loops and maintain their color with a certain probability or whether the probability that a node is part of a seed set is uniform versus non-uniform, independent versus correlated. For simplicity, we suppose here that nodes do not have self-loops and that the probability that a node is part of the seed set is indepedent and equal to $1/p_{\text{init}}$ for $0 < p_{\text{init}} < 1$.

\subsection{Formulation}

We are going to graph reconstruction problem from voter-model type cascades as a sparse recovery problem. We introduce the following notation:

\paragraph{Notation}
Denote by $(\vec m^t_k)_{1 \leq k \leq m; 1 \leq t \leq T}$ the set of blue nodes in cascade $k$ at time step $t$. In other words, $\vec m^t_k$ is a vector in $\mathbb{R}^n$ such that:

\[ (\vec m^t_k)_j = \left\{ 
  \begin{array}{l l}
    1 & \quad \text{if $j$ is \textsc{blue} in cascade $k$ at time step $t$}\\
    0 & \quad \text{otherwise}
  \end{array} \right.\]

Let $\vec x_i$ be the vector representing the neighbors of a node $i$ in the graph ${\cal G}$. In other words, $\vec x_i$ is a vector in $\mathbb{R}^n$ such that:

\[ (\vec x_{i})_j = \left\{ 
  \begin{array}{l l}
    1/\text{deg(i)} & \quad \text{if $j$ is a parent of $i$ in ${\cal G}$}\\
    0 & \quad \text{otherwise}
  \end{array} \right.\]

Let $w^t_{k,i}$ be the scalar representing the probability that in cascade $k$ at time step $t$ node $i$ stays or becomes \textsc{blue}.

\[ w^t_{k,i} = \left\{ 
  \begin{array}{l l}
    \frac{\text{Number of \textsc{blue} neighbors of i in $(k, t-1)$}}{\text{deg}(i)} & \quad \text{if $t\neq 0$}\\
    1/p_{\text{init}} & \quad \text{otherwise}
  \end{array} \right.\]

Our objective is to recover the neighbors for each node $i$ in ${\cal G}$. This corresponds to recovering the vectors $(\vec x_i)_{i \in {\cal G}}$. Notice that for each cascade $k$ and every time step $t < T - 1$, we have:

\begin{equation}
\langle \vec m_k^t, \vec x_i \rangle = \mathbb{P}(\text{node $i$ is \textsc{blue} in cascade k at $t+1$}) = w_{k,i}^{t+1} \nonumber
\end{equation}

We can concatenate each equality in matrix form $M$, such that the rows of $M$ are the row vectors $(\vec m^t_k)^T$, and for each node $i$, the entries of $\vec w_i$ are the corresponding probabilities:

$$M \vec x_i = \vec w_i$$

Note that if $M$ had full column rank, then we could recover $\vec x_i$ from $\vec b$. This is however an unreasonable assumption, even after having observed many cascades. We must explore which assumptions on $M$ allow us to recover $\vec x_i$ when the number of observed cascades is small. 

Further note that we do not immediately observe $\vec w_i$. Instead, we observe a bernoulli vector, whose j$^{th}$ entry is equal to $1$ with probability $(\vec w_i)_j$ and $0$ otherwise, corresponding to whether at the next time step $i$ became or stayed \textsc{blue}. We will denote this vector $\vec b_i$, and denote by $\vec \epsilon_i$ the zero-mean noise vector such that:

$$\vec w_i = \vec b_i + \vec \epsilon_i$$ 

For each node $i$, we can therefore reformulate the problem with our observables:

\begin{equation}
\label{eq:sparse_voter}
M \vec x_i = \vec b_i + \epsilon_i
\end{equation}

If the vector $\vec x_i$ is sufficiently sparse, i.e. node $i$ has sufficiently few neighbors, then we hope to apply common sparse recovery techniques to solve for $\vec x_i$.

\subsection{Example}


\begin{figure}
\centering
\includegraphics[width=0.3\textwidth]{images/voter.png}
\caption{A cascade in the voter model with time steps $t = 0,1,2$  over a graph with 5 vertices}
\end{figure}




To recover the neighbors of $v_1$, we get the following matrix $M$ for the example in Figure 1:

\begin{align*}
&\hspace{0.35cm} v_2 \hspace{0.2 cm}  v_3 \hspace{0.2 cm} v_4 \hspace{0.2 cm} v_5 &\\
\vspace{1 cm}
M  = &  \left( \begin{array}{cccc}
0 & 0 & 1 & 1 \\
1 & 1 & 0 & 0 \\
\end{array} \right)  &  \begin{array}{l} \hspace{ -4.5 cm} 
\text{time step 0} \\
 \hspace{ - 4.5cm}  \text{time step 1} \\
\end{array}
\end{align*}

and the vector $b_1$:

\begin{align*}
b_1  = &  \left( \begin{array}{c}
1 \\
1 \\
\end{array} \right) &  \begin{array}{l} \hspace{ - 5cm} 
\text{time step 1} \\
 \hspace{ - 5cm}  \text{time step 2} \\
 \end{array}
\end{align*}

\section{The Independent Cascade Model}
\label{sec:icc}
\subsection{Description}

In the independent cascade model, nodes have three possible states: susceptible, infected or active\footnote{The words `infected' and `active' are used interchangeably}, and inactive. All nodes start either as susceptible or infected. The infected nodes at $t=0$ constitute the seed set. Similarly to the voter model, we will suppose that nodes are part of the seed set with independent probability $p_{\text{init}}$. At each time step $t$, for every pair of nodes $(j,i)$ where $j$ is infected and $i$ is susceptible, $j$ attempts to infect $i$ with probability $p_{j,i}$. If $j$ succeeds, $i$ will become active at time step $t+1$. Regardless of $j$'s success, node $j$ will be inactive at time $t+1$: nodes stay active for only one time step. The cascade continues until time step $T_k$ when no active nodes remain.

\subsection{Formulation}

We are going to graph reconstruction problem from independent-cascade-model-type cascades as a sparse recovery problem. We introduce the following notation:

\paragraph{Notation}
We adapt slightly the definition of $(\vec m^t_k)$:

\[ (\vec m^t_k)_j = \left\{ 
  \begin{array}{l l}
    1 & \quad \text{if $j$ is \textsc{active} in cascade $k$ at time step $t$}\\
    0 & \quad \text{otherwise}
  \end{array} \right.\]

Similarly to the voter model, if $\vec p_i$ is the vector representing the `parents' of node $i$ in the graph ${\cal G}$ such that $(\vec p_i)_j$ is the probability that node $j$ infects $i$, then $\vec \theta_i$ is the vector:

\[ (\vec \theta_{i})_j = \left\{ 
  \begin{array}{l l}
    \log ( 1- p_{j,i}) & \quad \text{if $j$ is a parent of $i$ in ${\cal G}$}\\
    0 & \quad \text{otherwise}
  \end{array} \right.\]

Note that if $p_{j,i} = 0$, then $\theta_{j,i} = 0$ and we can write $\vec \theta = \log (1 - \vec p)$ element-wise. For $t>0$, let $w^t_{k,i}$ be the scalar representing the probability that in cascade $k$ at time step $t$ node $i$ becomes \textsc{active}, conditioned on the state of cascade  $k$ at time step $t-1$. As long as node $i$ has not yet been infected in cascade $k$, we can write as in the voter model, we can observe that: 

$$\langle \vec m_k^t, \vec \theta_i \rangle = \vec w^t_{k,i}$$

As in the voter model, by concatenating all row vectors $\vec m^t_{k,i}$ such that $t < t^k_i$, where $t^k_i$ is the infection of node $i$ in cascade $k$, in the matrix $M_i$, we can write:

$$M_i \vec \theta_i = \vec w_i$$

Once again, we do not have access to $\vec w_i$. Instead, for $t<T_k$, we have access to the bernoulli variable $b_i$ of parameter $1 - e^{w_i}$, indicating whether or not node $i$ became active at time step $t+1$ in cascade $k$. Recovering $\vec \theta_i$ for every node $i$ can be written as the following problem:

\begin{equation}
e^{M_i \vec \theta_i} = 1 - \vec b_i + \vec \epsilon_i
\end{equation}

Note that this is not exactly a sparse recovery model due to the non-linear exponential term. It is of the author's opinion however that if the probabilities $p_{j,i}$ are restrained to a bounded interval $[0, 1- \eta]$, then most of the results which hold for the linear voter model will continue to hold in this case.

\subsection{Example}

\begin{figure}
\centering
\includegraphics[width=0.3\textwidth]{images/icc.png}
\caption{A cascade in the independent cascade model with time steps $t = 0,1,2$  over a graph with 5 vertices}
\end{figure}



To recover the neighbors of $v_5$, we get the following matrix $M$ for the example in Figure 2:
\begin{align*}
&v_1 \hspace{0.2 cm}  v_2 \hspace{0.2 cm} v_3 \hspace{0.2 cm} v_4 \\
\vspace{1 cm}
M  = &  \left( \begin{array}{cccc}
1 & 0 & 0 & 0 \\
0 & 1 & 1 & 0 \\
\end{array} \right)   \begin{array}{l} \hspace{ 1cm} 
\text{time step 0} \\
 \hspace{ 1cm}  \text{time step 1} \\
 \end{array}
\end{align*}

and the vector $b_5$:

\begin{align*}
b_5  = &  \left( \begin{array}{c}
0 \\
1 \\
\end{array} \right)  \begin{array}{l} \hspace{ 1cm} 
\text{time step 1} \\
 \hspace{ 1cm}  \text{time step 2} \\
 \end{array}
\end{align*}



\section{Sparse Recovery}

\subsection{Introduction}
The literature on sparse recovery has expanded rapidly since a series of papers were published by Candès, Tao and Romberg in 2006. In 8 years, over 7000 articles have been written on the topic since these few seminal publications. Many of them show that if the matrix $M$ is `well-conditioned', then any $s$-sparse signal can be recovered from ${\cal O}(s \log n)$ measurements (a row in the matrix $M$). We explored several of these papers to understand which conditions on $M$ were most reasonable.

\subsection{Exploration}
\paragraph{Disclosure}
After an unsatisfactory attempt at adapting \cite{candestao}, we cite the main results of \cite{candestao} and give arguments as to why we believe, with more work, these assumptions can be made reasonable in the case our model. Due to the limited time frame of this project, we were unable to find a complete and rigorous argumentation as to why these assumptions hold. The following paragraphs should be understood as motivation for this line of research rather than concrete guarantees.

\paragraph{Framework}
Suppose that each row $m$ of $M$ is taken from a distribution ${\cal F}$, such that $\mathbb{E}(m m^T) = I_n$. Suppose further that $\|m\|_{\infty} \leq \mu$ and that $\|M^T \epsilon\|_{\infty} < 2.5 \sqrt{\log n}$. For every node $i$, consider the following optimization program:

\begin{equation}
\label{eq:optimization_program}
\min_{x_i \in \mathbb{R}^n} \frac{1}{2}\|M \vec x_i - b_i\|^2_2 + \lambda \| x_i\|_1
\end{equation}

\begin{theorem}
\label{thm:candes}
If the number of rows $l$ of our matrix $M$ verifies $l \geq C \mu s \log n$ (where $C$ and $s$ are constants), then with high probability and for the right choice of $\lambda$, the solution $\hat x$ to Eq.~\ref{eq:optimization_program} obeys:
\begin{equation}
\label{eq:first_guarantee}
\| \hat x - x^* \|_2 \leq C (1 + \log^{3/2}n)\left[\frac{\|x^* - \tilde x\|_1}{\sqrt{s}} + \sigma \sqrt{\frac{s\log n}{l}}\right]
\end{equation}
where $\sigma$ is a constant dependent on the variance of $\epsilon$, and $\tilde x$ is the best s-sparse approximation to the ground truth $x^*$. In other words if node $i$ is of degree at most $s$, then $\|x^* - \tilde x\|_1 = 0$ since $x^*$ is $s-sparse$.\footnote{The result in the paper supposes that $x_i^*$ is $s$-sparse, i.e. that node $i$ has at most $s$ neighbors.}
\end{theorem}

\paragraph{Discussion}
\begin{itemize}
\item Note how these results only hold for the voter model. For the independent cascade, our first unsucessful attempt at applying \cite{candestao} showed that its results could be extended to the non-linear case with different constants under the assumption that $\forall i,j, \ p_{j,i} < 1 - \eta$. Since the spirit of the proofs are very similar between the two papers, we hope that Theorem~\ref{thm:candes} can also be extended.
\item The assumption $\|M^T \epsilon\|_{\infty} < 2.5 \sqrt{\log n}$ can most likely be shown to hold in our framework. Indeed, regardless of the normalization of $M$, $\epsilon$ will be of mean 0. Therefore, by a Chernoff-like bound followed by a union bound, $\|M^T \epsilon\|_{\infty}$ will be small with high probability.
\item The more difficult assumption to adapt is $\mathbb{E}(m m^T) = I_n$. The authors mention that their results can be extended to cases where $\mathbb{E}(m m^T) \approx I_n$. Suppose that restrict ourselves to the first row of every cascade. In that case, the entries of our matrix are independent bernoulli variables equal to $1$ with probability $p_{\text{init}}$ and 0 otherwise. If we normalize our matrix with $1/\sqrt{p_\text{init}}$, then $\mathbb{E}(m' m'^T)$ is an $n \times n$ matrix with 1 on the diagonal on $p_{\text{init}}$ everywhere else, which is $\approx I_n$ if $p_{\text{init}}$ is small.
\item The authors also mention that if the restricted isometry property (RIP) holds with $\delta < \frac{1}{4}$, then a similar theorem could easily be shown. The restriced isometry property characterizes a quasi-orthonormality of the measurement matrix M on sparse vectors: For all k, we define $\delta_k$ as the best possible constant such that for all unit-normed ($l$2) and k-sparse vectors $x$:

\begin{equation}
1-\delta_k \leq \|Mx\|^2_2 \leq 1 + \delta_k \nonumber
\end{equation}

In general, the smaller $\delta_k$ is, the better we can recover $k$-sparse vectors $x$. In the next section, we explore the $\delta$ constant of such a matrix, and show that in certain contexts it can be shown to be $\delta < \frac{1}{4}$.
\end{itemize}

\section{Experiments}
\label{sec:exp}

\subsection{Verifying the RIP constants}

\paragraph{Motivation} Though it is not yet clear which normalization to choose for the columns of matrix $M$, nor whether we should include only the first step of every cascade or the cascade in its entirety, we ran simple experiments to show that the measurement matrix $M$ has `well-behaved' RIP constants. Indeed, the authors of \cite{candes} mention in their proofs that with an RIP strictly less than $\frac{1}{4}$, the assumptions of the theorem would hold, with slightly different constants. Due to a lack of time, we were unable to verify these claims or understand what the final statement of the theorem would be. In general, a measurement matrix with small RIP constants is known to recover sparse vectors well. These experiments should be considered as another reason to think that the sparse recovery framework is appropriate to the graph reconstruction problem.  

\paragraph{In practice} Finding a non-trivial RIP constant is strongly NP-hard, and as such we were unable to test on networks with more than 25 nodes. To compute the $k-th$ RIP constant ($\delta_k$ with our notation), we must extract $k$ columns of $M$, and compute the following programs on the resulting matrices $(M_k)$ (there are $\binom{n}{k}$ of them!):

\begin{table}
\centering
\begin{tabular}{c | c | c | c }
c = 100 & c = 1000 & c = 100 & c = 1000 \\
$p_\text{init}$ = .1 & $p_\text{init}$ = .1 & $p_\text{init}$ = .05 & $p_\text{init}$ = .05 \\
\hline
$\delta_4 = .54$ & $\delta_4 = .37$ & $\delta_4 = .43$ & $\delta_4 = .23$ \\
\end{tabular}
\caption{Results on dataset of 22 nodes extracted from the Karate club network. $p_\text{init}$ is the probability that a node is part of the seed set. $c$ is the number of cascades. $\delta_4$ is the RIP constant for $4$-sparse vectors. Cascades are generated using the independent cascade model.}
\end{table}


\begin{align}
\label{eq:lower_bound}
\min \qquad & \| M_k \vec x \|_2 \\
\| \vec x\|_2 = 1 & \nonumber
\end{align}
\begin{align}
\label{eq:upper_bound}
\max \qquad & \| M_k \vec x \|_2 \\
\| \vec x\|_2 = 1 & \nonumber
\end{align}

\paragraph{Discussion} Note how choosing the minimum of Eq~\ref{eq:lower_bound} over all extracted matrices returns the lowest possible $1 - \delta^-_k$ and choosing the maximum of Eq~\ref{eq:upper_bound} over all extracted matrices returns the highest possible $1 + \delta^+_k$. The smallest admissible $k-$RIP constant is therefore: $\max(\delta^+_k, \delta^-_k)$. Equations \ref{eq:lower_bound} and \ref{eq:upper_bound} can be solved efficiently since they correspond to the smallest and greatest eigenvalues of $M^TM$ respectively. 

The results of our findings on a very small social network (a subset of the famous Karate club), show that as the number of cascades increase the RIP constants decrease and that if $p_\text{init}$ is small then the RIP constant decrease as well. Finally the constants we obtain are either under or close to the $.25$ mark set by the authors of \cite{candes}. By choosing appropriate parameters (number of cascades large enough, and $p_\text{init}$ small enough), we can place ourselves in the framework of this theorem for nodes of degree $4$ in a dataset of size 22.


\subsection{Testing our algorithm}

\begin{figure}
\centering
\label{fig:comparison-with-greedy}
\includegraphics[scale=.35]{images/greedy_sparse_comparison.pdf}
\caption{Plot of the F1 score for different number of cascades for both the \textsc{greedy} algorithm and Algorithm~\ref{eq:optimization_program}. The cascades, graphs, and edge probabilities were identical for both algorithms. The dataset is a subset of 333 nodes and 5039 edges taken from the Facebook graph, made available by \cite{snap}. We chose $p_\text{init}=.05$, and the edge probability were chosen uniformly between $0$ and $0.8$. For Algorithm~\ref{eq:optimization_program}, we chose $\lambda=10$ and kept all edges with probability greater than $.1$ as true edges.}
\end{figure}

We were only able to test our algorithm in the case of the independent cascade model and compare it to the \textsc{greedy} baseline suggested in \cite{netrapalli}. We picked a small subset of 333 nodes and 5039 edges from the Facebook graph, made available by \cite{snap}. As we can see from Figure~\ref{fig:comparison-with-greedy}, our algorithm outperforms considerably the \textsc{greedy} baseline. Of things to take into consideration are the fact that the graph is very small, that the \textsc{greedy} algorithm has no guarantee for non-tree networks\footnote{Their authors do mention that it performs well in practice on non-tree networks as well}, and that \textsc{greedy} is a very simple heuristic. To get an accurate picture of our algorithm's contribution, it would be necessary to compare it to every other algorithm suggested in \cite{netrapalli} and \cite{abrahao}.


\section{Conclusion}

In conclusion, we described how a well-studied framework could be applied to the graph reconstruction problem. We give arguments as to why applying this framework could lead to better guarantees than pre-existing results. We also test our model in practice and show that it performs better than a simple heuristic, validating our approach to a certain extent.

\begin{thebibliography}{42}

\bibitem{gomez}
Gomez-Rodriguez, M., Leskovec, J., and Krause, A.
\newblock {\it Inferring Networks of Diffusion and Influence}
\newblock In Proc. 16th ACM SIGKDD international conference on Knowledge discovery and data mining, KDD'10, pages 1019--1028,
\newblock 2010

\bibitem{netrapalli}
Netrapalli, P., and Sanghavi, S.,
\newblock {\it Learning the Graph of Epidemic Cascades}
\newblock SIGMETRICS Perform. Eval. Rev., 40(1): 211--222,
\newblock 2012


\bibitem{abrahao}
Abrahao, B., Chierichetti, F., Kleinberg, R., and Panconesi, A.
\newblock{\it Trace Complexity of Network Inference}
\newblock In Proc. 19th ACM SIGKDD international conference on Knowledge discovery and data mingin, KDD'13, pages 491--499,
\newblock 2013

\bibitem{candestao}
Candès, E. J., Romberg, J. K., and Tao, T.
\newblock {\it Stable Signal Recovery from Incomplete and Inaccurate Measurements},
\newblock Comm. Pure Appl. Math., 59: 1207--1223,
\newblock 2006,

\bibitem{candes}
Candès, E., and Plan, Y.
\newblock {\it A Probabilistic and RIPless Theory of Compressed Sensing}
\newblock Information Theory, IEEE Transactions on, 57(11): 7235--7254,
\newblock 2011.

\bibitem{snap}
Leskovec, J., and Krevl, A.,
\newblock {\it SNAP Datasets: Stanford Large Network Dataset Collection},
\newblock 2014.


\end{thebibliography}

\end{document}