Probability of global epidemic Consider $\mathbb{Z}^2$ as a graph, where each node has four neighbours. 4 signals are emitted from $(0,0)$ in each of four directions (1 per direction) . A node that receives one signal (or more) at a timestep will re-emit it along the 4 edges to its four neighbours at the next time step. A node that did not receive a signal at the previous timestep will not emit a signal irrespective of whether it earlier received a signal. There is a $50\%$ chance that a signal will be lost when travelling along a single edge between two neighbouring nodes.
A node that receives more then 1 signal acts the same as if it received only 1.
The emitting of a signal in each of the 4 directions are independent events.
What is the probability that the signal will sometime arrive at $(10^5,10^5)$?
Research: Simulations show: Yes. ~90%
What is the probability that the signal will sometime arrive at $(x,y)$ if a signal traveling along an edge dies with probability $0<p<1$?
What is the least p, for which the probability that N initial random live cells die out approaches 0, as N approaches infinity? Experiment shows p close to 0.2872.
In $\mathbb{Z}^1$, $p_{min}=0.6445...$, how to calculate this?
In $\mathbb{Z}^3$, $p_{min}=0.1775...$, how to calculate this?
 A: This does not compute an exact value but some rigorous upper bounds of the probability $r$ that every site in $\mathbb Z^2$ is reached, eventually. 
Consider the Galton-Watson branching process where the progeny of every individual has the distribution of the number of signals transmitted at time $n+1$ by any site reached at time $n$. Thus the progeny of any node in the Galton-Watson tree is $0$, $1$, $2$, $3$ and $4$ with respective probabilities $1$, $4$, $6$, $4$ and $1$ divided by $16$. 
Then the generation $n$ of the process overestimates the number of sites reached at time $n$ since the signals received from different neighbors stay separated. Hence the probability $q$ that the tree becomes extinct underestimates the probability that the cloud of sites receiving a signal becomes empty. On the other hand, if the cloud survives forever, the path of any eternally transmitted signal is a simple random walk on $\mathbb Z^2$, and these are recurrent hence every site is reached by the trace of the cloud on the grid, eventually.
The extinction probability $q$ is classically the smallest root in $[0,1]$ of the one-step equation 
$$
16q=1+4q+6q^2+4q^3+q^4,
$$ 
which is $q=0.087378$. Hence every site on the grid is reached by the original transmission process with probability $r\lt91.2622\%$.
Now, the probability $r_T$ that $T=(10^5,10^5)$ is reached corresponds to the fact that generation $n=2\cdot10^5$ of the tree is not empty, whose probability $1-q_n$ is very close to $1-q$, hence $r_T\lt1-q_n\approx1-q$.
Nota: One gets a simpler bound noticing that, if the initial site transmits nothing, no other site is reached, and this happens with probability $1/16=6.25\%$, hence $r\lt93.75\%$.
A: In the same vein as Didier's answer providing bounds on the extinction probability $q$, we can also obtain bounds on the transmission probability $p$ required for the probability of $N$ signals dying out to go to $0$ as $N\to\infty$, which is the probability required for the extinction probability $q$ not to be $1$.
Didier's equation for $q$ can be rewritten as $q=(1-(1-q)/2)^4$, which says that the signal goes extinct if all four signals die, either because they don't make it across the edge ($1/2$) or because they go extinct afterwards ($1-q$). Generalizing to $d$ dimensions and transmission probability $p$, this is
$$q=(1-p(1-q))^{2d}\;.$$
This equation always has a root at $1$. For sufficiently large $q$, it also has a second root in $[0,1]$, and the extinction probability is given by that root. The critical case in which the extinction probability becomes $1$ occurs when these two roots coincide. Differentiating with respect to $q$ and substituting $q=1$ yields the condition for $1$ to be a double root:
$$
\begin{eqnarray}
1&=&2dp(1-p(1-q))\;,
\\
1&=&2dp\;,
\\
p&=&\frac1{2d}\;.
\end{eqnarray}
$$
This is a lower bound, since in the case of confluent signals there are fewer signals to keep the fire burning. Your experiments seem to indicate that the bound becomes better with increasing dimension, which makes sense since the signals become less likely to collide.
A: This is not an answer, but was too long for a comment.
More like some observations and additional, related questions
on the theme of whether generating functions would help.
My particular choice of generating functions model additive signals
and therefore could be used as upper bounds for nonadditive ones.
My approach below is probably not original.
Good combinatorial interpretations probably already exist.
This is a Galton-Watson stochastic branching process.
The one-dimensional case has generating function
in the formal power series ring with indeterminates $x,x^{-1}$
(or its quotient with the ideal $\left< x^{-1}x-1 \right>$)
given by
$$
\eqalign{
p_t^*(x)
&= p^t \left(x+x^{-1}\right)^t
 = p^t \sum_{s=0}^{t} {t \choose s} x^{2s-t}
 = p^t \sum_{s=0}^{t} {t \choose s} x^{t-2s}
\cr &
= p^t
  \left(
    \left\lfloor\tfrac{t}2\right\rfloor-
    \left\lfloor\tfrac{t-1}2\right\rfloor
  \right)
  {t \choose t/2}
+ p^t
  \sum_{s=0}^{\left\lfloor\tfrac{t-1}2\right\rfloor}
  {t \choose s}
  \left(x^{t-2s}+x^{2s-t}\right)
}
$$
for time $t$, i.e. the coefficient of $x^m$
is the probability the signal reaches coordinate
$x=m$ at time $t$. The first term on the second line,
the constant term, representing the probability
of revisiting the origin, is positive with
maximum coefficient for $t$ even, but $0$ for $t$ odd.
The differences in parentheses above,
$$
e(t)=
 \left\lfloor\tfrac{t}2\right\rfloor-
 \left\lfloor\tfrac{t-1}2\right\rfloor
=
 \tfrac{t+1}2-2\left\lfloor
 \tfrac{t+1}2\right\rfloor
=\frac{1+(-1)^{t+1}}{2}
=(t+1)~\text{mod}~2

$$
indicates evenness.
We could also take the last expression and throw away
the (strictly) negative powered terms $x^{2s-t}$,
redefining the coefficient of $x^m$ to be the probability
that a signal starting at the origin at time $0$
reaches $x=\pm m$ (one or the other) after time $t$,
and call this the reduced generating function
$$
p_t(x)
= p^t \left[ e(t) {t \choose t/2}
+ \sum_{s=0}^{\left\lfloor\tfrac{t-1}2\right\rfloor}
  {t \choose s} x^{t-2s}
\right].
$$
Note in particular that $e(t){t\choose t/2}$ is therefore the
number of paths of length $t$ which start and end at the origin.
Summing over $t\ge0$, we get a generating function
$$
p(x)=\sum_{t=0}^{\infty} p_t(x)
$$
for the expected amount of time spent at each coordinate,
i.e. the coefficient of $x^m$ is the sum of probabilities
that the signal is at one of $x=\pm m$ over all times $t\ge 0$.
Of particular interest to us is
the expected total time spent at the origin,
$$
p(0)=\sum_{t=0}^\infty{2t\choose t}p^{2t}.
$$
Combinatorial note: there is no double-counting
(no need for P.I.E.) because each term counts only
the terminus point over all paths of length $t$.
This acts like a "driving force" of a
harmonic oscillator damped by a probability of dying out.
By Stirling's approximation and the root test,
this converges, for $p<\frac12$.

Claim $\left(\mathbb{Z}^1\right)$:
  For $r=2p<1$, the expected total number of
  visits to the origin converges to
  $$p(0)=\frac{1}{\sqrt{1-r^2}}.$$

Proof: this follows from the general binomial theorem,
as in e.g. Yuan 2001 proposition 4.3:
$$
\left(-4\right)^n {-\tfrac12 \choose n} =
\left(-4\right)^n \frac{(-1)(-3)\cdots(1-2n)}{2^n\,n!} =
\frac{(2n)!}{(n!)^2}.
$$
Now for the two-dimensional case, we can avoid using the
multinomial theorem expansion, if we like, with
$$
\eqalign{
p_t^*(x,y)
&= p^t \left(x+x^{-1}+y+y^{-1}\right)^t
= \sum_{s=0}^{t} {t \choose s} p_s^*(x) \, p_{t-s}^*(y)
\cr
&= p^t
  \sum_{s=0}^{ t }
  \sum_{u=0}^{ s }
  \sum_{v=0}^{t-s}
  { t \choose s}
  { s \choose u} 
  {t-s\choose v}
  x^{2u-s} y^{2v+s-t}
}
$$
or, throwing away the negative powered terms,
$$
\eqalign{
p_t(x,y)
&= \sum_{s=0}^{t} {t \choose s} p_s(x) \, p_{t-s}(y)
%%
%% This formula could be simplified analogously to p_t^*(x) above:
%%
%\cr
%&= p^t
%  \sum_{s=0}^{ t }
%  \sum_{u=0}^{ s }
%  \sum_{v=0}^{t-s}
%  { t \choose s}
%  { s \choose u} 
%  {t-s\choose v}
%  x^{2u-s} y^{2v+s-t}
}
$$
and take the analogous sum $p(x,y)$.
Note that $p_t(0,0)=0$ for odd $t$, and that (for $2t$ even)
$$
p_{2t}(0,0)
=p^{2t}\sum_{s=0}^t{2t\choose 2s}{2s\choose s}{2t-2s\choose t-s}
=p^{2t}\sum_{s=0}^t{t\choose s}^2
=p^{2t}{2t\choose 2}^2
$$ 
which is $p^t$ times the square of a 
central binomial coefficient; cf. sequence A000984.
(The latter two formulas, with the squared binimial
coefficients, are offered here without proof.)
Hence
$$
p(0,0)=\sum_{t=0}^\infty{2t\choose t}^2p^{2t},
$$
which will converge for $p<\frac14$.

Claim $\left(\mathbb{Z}^2\right)$:
  For $r=4p<1$, the expected total number of
  visits to the origin converges to
  $$p(0,0)={}_2F_1\left(\tfrac12;\tfrac12;1;r^2\right)=\frac1{M(1,\sqrt{1-r^2})}$$
  converges to a special case
  of the Gaussian hypergeometric series ${}_2F_1$
  which is a complete elliptic integral of the first kind and
  can be calculated with extremely rapid convergence
  (Gauss 1799)
  by the arithmetic-geometric mean
  $$M(a,b)=\int_0^{\frac\pi2}\frac{d\phi}{\sqrt{(a\cos\phi)^2+(b\sin\phi)^2}}.$$

The first OP question is perhaps related to the expected time
$$
c_{mn} =
\mathbb{E} \left[ T_{mn} \right] =
\frac1{m!\,n!}
\left.
\left(
\partial_x^{\,m}
\partial_y^{\,n}
 p(x,y)
\right)
\right|_{(0,0)}
\qquad(m,n\ge0)
$$
spent at one of the four points $(\pm m,\pm n)$
which is the sum of the $x^my^n$ terms over all 
$p_t(x,y)$ and can be computed deterministically
with no further simplification in $O(m^3n^3)$ time
by adding products of binomial coefficients.
Perhaps there are ways to simplify this
using some combination of binomial identities,
parity and symmetry, or with asymptotics,
to get a reasonable approximation
in a feasible amount of time.
The familiar formula with the higher order derivative
using the Taylor series trick for the desired coefficient
is valid on the reduced generating function for $m,n\ge 0$.
As a wild speculation, is the answer
to the (first and) second OP question
then
$$
\frac{c_{mn}}{1-p(0,0)}
$$
or something in a similar vein, analogously to the
gambler's ruin problem or the reasoning of @joriki or @Didier?
Conjectures for higher dimensions:

Conjecture: The expected number of visits to the origin
  in $\mathbb{Z}^d$ is
  $$ \sum_{n=0}^\infty {2n \choose n}^d p^{2n}. $$
Conjecture: The number of paths in $\mathbb{Z}^d$ of length $2n$
  which start and end at the origin is
  $$ {2n \choose n}^d. $$
  Note: $d=2$ is Richard Stanley's
  bijective proof problem #230. 

A: This may be too simplistic, and I only get a $71.0\%$ result.  I assume that around the node $(10^5,10^5)$ all the neighboring nodes have the same probability of being reached, call it $P$ to distinguish it from the $p$ in the question.  
One other assumption:  transmissions works $1/2$ the time.  
        *

   *    *    *
       (a,a)

        *

So $(a,a)$ could receive it from any of the four neighbors.  
At each of the neighbor nodes I use the conditional probability that the node receives the signal before the center node $(a,a)$ does.  I assume that  equals $.75P$ since the neighbor node then did not receive the signal from $(a,a)$.
Then I have the following equation using the binomial probability formula
$$\binom{4}{1} (.75 P) (1-.75P)^3(1-0.5)+ \binom{4}{2} (.75P)^2 (1-.75P)^2 (1-.25)+ \binom{4}{3} (.75)^3 (1-.75P) (1-.125)+ \binom{4}{3} (.75)^4 (1-.0625) = P$$
Solving, $P \approx .703$.
This approach fails when close to the origin, I think.
