Seed Selection in the Heterogeneous Moran Process
Petros Petsinis Andreas Pavlogiannis Josef Tkadlec Panagiotis Karras


Symmetries in particle physics:
 from nuclear isospin to the quark model
    Giovanna F. do Valle
    April 29, 2024
========================================================================






The Moran process is a classic stochastic process that models the rise and takeover of novel traits in network-structured populations.
In biological terms, a set of mutants, each with 
fitness ∈(0,∞) invade a population of residents with fitness 1. 
Each agent reproduces at a rate proportional to its fitness and each offspring replaces a random network neighbor. 
The process ends when the mutants either fixate (take over the whole population) or go extinct.
The fixation probability measures the success of the invasion.
To account for environmental heterogeneity,
we study
a generalization of the Standard process, called the Heterogeneous Moran process.
Here, the fitness of each agent is determined both by its type (resident/mutant) and the node it occupies.
We study the natural optimization problem of seed selection: given a budget k, which k agents should initiate the mutant invasion to maximize the fixation probability?
We show that the problem is strongly inapproximable: it is -hard to distinguish between maximum fixation probability 0 and 1.
We then focus on mutant-biased networks, where each node exhibits at least as large mutant fitness as resident fitness.
We show that the problem remains -hard, but the fixation probability becomes submodular, and thus the optimization problem admits a greedy (1-1e)-approximation.
An experimental evaluation of the greedy algorithm along with various heuristics on real-world data sets corroborates our results.






§ INTRODUCTION


Modeling and analyzing the spread of a novel trait (e.g., a trend, meme, opinion, genetic mutation) in a population is vital to our understanding of many real-world phenomena. 
Typically, this modeling involves a network invasion process: nodes represent agents/spatial locations, edges represent communication/interaction between agents,  and local stochastic rules define the dynamics of trait spread from an agent to its neighbors.

Network diffusion processes raise several optimization challenges, whereby we control elements of the process to achieve a desirable emergent effect. 
A well-studied problem is that of influence maximization, which calls to find a seed set of agents initiating a peer-to-peer influence dissemination that maximizes the expected spread thereof; the problem arises in various diffusion models, such as Independent Cascade and Linear Threshold <cit.>, the Voter model <cit.>,  and geodemographic models of agent mobility <cit.>.


Diffusion processes also play a key role in evolutionary dynamics, which model the rules underpinning the sweep of novel genetic mutations in populations and the emergence of new phenotypes in ecological environments <cit.>. 
A classic evolutionary process is the Moran process <cit.>. 
In high level, a set of mutants, each  with fitness ∈ (0,∞), invade a preexisting population of residents, each with fitness normalized to 1. 
Over time, each agent reproduces with rate proportional to its fitness, while the produced offspring replaces a random neighbor. 
In the long run, the new mutation either fixates in the population (i.e., all agents become mutants) or goes extinct (i.e., all agents remain residents). 
The probability of fixation is the main quantity of interest, especially under advantageous mutations (>1).

Network structure affects the fixation probability <cit.>, and may both amplify it <cit.> and suppress it <cit.>, while certain structures nearly guarantee mutant fixation <cit.>.
The Moran process thus provides a simple stochastic model by which a community of communicating agents reaches consensus; one option has an advantage over another, yet its prevalence (i.e., fixation) depends on the positioning of its initial adherents (i.e., mutants) and on the network structure.


Recent work aims to make the Moran process more realistic by incorporating some form of environmental heterogeneity <cit.>.
Here, the fitness of an agent is not only a function of its type (resident/mutant), but also of its location in space, i.e., the node that it occupies.
For example, in a biological setting, the ability to metabolise a certain sugar boosts growth more in environments where such sugar is abundant.
Similarly, in a social setting, the spread  a trait is more, or less viral depending on the local context (e.g., ads, societal predispositions).
Analogous extensions have been recently considered for the Voter evolutionary model <cit.>.



In this work we generalize the classic Moran process to account for complete environmental heterogeneity, obtaining the Heterogeneous Moran process: for every network node u, a mutant (resp., resident) occupying u exhibits fitness (u) (resp., (u)) specific to that node.
We then study the  natural optimization problem of seed selection: given a budget k, which k nodes should initiate the mutant invasion so as to maximize the fixation probability?
Although the seed selection problem has been studied extensively in other diffusion models,
this is the first paper to consider it in Moran models.
We obtain upper and lower bounds for the complexity of this problem in our Heterogeneous model,
which also imply analogous results to other relevant Moran models.



Contributions
Our main theoretical results are as follows (see <ref> for a summary in the context of Moran models).


  * We prove that computing the fixation probability admits a FPRAS on undirected and unweighted networks that are mutant-biased, where (u)≥(u) for every node u. 

  * We show that the optimization problem is strongly inapproximable: for any 0<ε<12, it is -hard to distinguish between maximum fixation probability ≤ε and >1-ε.

  * We then focus on mutant-biased networks.
We show that the optimization problem remains -hard to solve exactly, but the fixation probability becomes submodular, yielding
a greedy (1-1e)-approximation.

Further, we evaluate the greedy algorithm and some standard heuristics for seed selection on real-world data.
Our experiments indicate that the
greedy algorithm outperforms all heuristics and uncovers 
high-quality seed sets for varying real-world datasets and problem parameters.
Due to space constraints, some proofs appear in the <ref>.

Technical Challenges
The problem of seed selection was studied recently under the Voter model, which bares some resemblance to the Moran model <cit.>.
However, the two models are distinct, and results in one do not transfer to the other.
Some novel technical challenges we address in this work are as follows.


  * Our NP-hardness and inapproximability proofs are fundamentally different from the NP-hardness of <cit.>, and are not limited to weak selection (mutant advantage ϵ 0).

  * Our submodularity proof is based on introducing a novel variant of the Moran process that we call the Loopy process.
This also allows us to show that the Heterogeneous Moran process
is a special case of the Two-Graph Moran process <cit.>, thereby extending
our hardness results to the latter.

  * Our model accounts for environmental heterogeneity,
while the Voter model in <cit.> does not. 
This complicates our proof for FPRAS.






§ PRELIMINARIES

In this section we introduce the Heterogeneous Moran process and the problem of seed selection.

Population structure 
We consider a population of individuals arranged in space. 
The population structure is represented as a weighted directed graph G =(V, E, ), 
where each node u∈ V represents a single agent, 
each edge (u,v)∈ E represents the fact that u interacts with (influences) v, 
and (u,· ) is a probability distribution capturing the frequency at which u influences v. 
We require that G is strongly connected, i.e., any two nodes are connected by a sequence of edges of non-zero weight.
We call G undirected if E is symmetric and (u,· ) is uniform.

Fitness graphs 
Trait diffusion in the Heterogeneous Moran process occurs by associating each node with a type: at each moment in time, each node u is either resident or mutant. 
Moreover, u is associated with a fitness that is type-dependent and represents the rate at which u influences its neighbors while being resident or mutant.
We denote the respective fitness values by (u) and (u); concretely, these are functions ,  V→ (0,∞).
We call the triplet =(G, (, )) a fitness graph, and denote the minimum and maximum resident and mutant fitnesses in  as
=min_u∈ V(u), and
=max_u∈ V(u).
We call  mutant-biased if for all u∈ V, we have (u)≥(u).
 
The Heterogeneous Moran process 
A configuration is a subset of nodes ⊆ V, representing the mutant nodes in  at some time point.
The fitness of node u in  is defined as 


    _(u) = 
    (u),   if  u∈
    (u),   otherwise


i.e., it is (u) if u is mutant and (u) if u is a resident.
At time t=0 a seed set ⊆ V 
specifies the nodes where mutant invasion begins.
The Heterogeneous Moran process is a 
discrete-time 
stochastic process _0, _1, …, of stochastic configurations _t ⊆ V, where
_0 = and 
for each t>0,
_t+1 is obtained from _t by two successive random steps:


  * Birth Event: Pick a node u for reproduction with probability proportional to its fitness,
_(u)/∑_v∈ V_(v).

  * Death Event: Pick a neighbor v of u with probability (u,v) and make v have the same type as u.

Note that the mutant set can both grow and shrink over time.

<ref> illustrates the process on a small example.



Relation to other Moran processes
We recover the Standard Moran process <cit.>
as a special case of the Heterogeneous process with (u)=1 and (u) is a constant for all u∈ V, expressing the mutant fitness.
The Neutral Moran process is a further special case of the Standard process, having δ=0 (i.e., residents and mutants have equal fitness).
The Positional Moran process <cit.>
is a parameterization of the Standard process by an active set of nodes ,
which define the node fitness as _(u)=1+δ if u∈∩ and _(u)=1 otherwise.
This is also a special case of the Heterogeneous process we study here, with (u)=1 and (u)=1+δ if u∈ and (u)=1 otherwise.
The Two-Graphs Moran process <cit.> extends the Standard model by accounting for heterogeneity in mobility, as opposed to fitness: mutants and residents propagate via different, type-specific graphs G_M and G_R, respectively, over the same set of nodes but with  different edge sets.
The Two-Graphs process is more general than our Heterogeneous process; though this connection is not obvious, it is formally implied by an intermediate result we establish in <ref> towards submodularity.

Fixation probability 
In the long run, the process reaches a consensus state _t = V (mutant fixation) or _t=∅ (mutant extinction). 
The fixation probability _() is the probability that fixation occurs when the process runs on a fitness graph  with seed set .
The complexity of computing the fixation probability is an open problem, even for the Standard process.
In the next section we prove that the fixation probability can be approximated efficiently via Monte Carlo simulations on mutant-biased, undirected fitness graphs.



The seed-selection problem 
The standard optimization question in invasion processes is optimal seed placement: given a budget k, which k nodes ^* should initiate the mutant invasion so as to maximize the fixation probability?


    ^*=_⊆ V, ||≤_().



The optimal seed depends on the graph structure,  budget , and node fitnesses.
<ref> showcases this intricate relationship, even when all residents have fitness 1.
In particular, the optimal seed ^* may contain
(i) only nodes with the largest mutant fitness 
(=3, left),
(ii) nodes with both large and small mutant fitness (=3, middle), or
(iii) even only nodes with the smallest mutant fitness (=3, right).
Moreover, the optimal seed is not monotonic on , i.e., increasing  may yield an optimal seed set that is not a superset of the previous one, and the two may even be disjoint (left, =1 vs =3).





§ COMPUTING THE FIXATION PROBABILITY

Given a seed set  and a fitness graph =(G,(,)), the complexity of computing the fixation probability _() is a long-standing open problem even in the Standard Moran process.
This is in sharp contrast to standard cascade models of influence spread, for which the spread function can be approximated efficiently <cit.>.

In the neutral setting ((u)=(u) for all u), 
the fixation probability is linear, _() = ∑_u∈_({u}). 
When the graph is also undirected, _() = ∑_u ∈1(u)/∑_v∈ V1(v), where (x) is the degree of  x <cit.>. 
On the other hand, no closed-form solution is known for the non-neutral setting.
However, on undirected graphs the expected time until convergence is polynomial, yielding a fully polynomial-time randomized approximation scheme (FPRAS) via Monte Carlo simulations <cit.>.
The following lemma generalizes the above result to the Heterogeneous process on mutant-biased graphs.
This is in sharp contrast to non-biased graphs, on which the expected time is exponential in general <cit.>.

lemmalemexpectedtime
Given an undirected and mutant-biased fitness graph  and a seed set ⊆ V,
the expected time to convergence T(,) satisfies T(,)≤(n^2·/)^3.


For a configuration , we define the potential function Φ(X)=∑_u∈(u)/(u), where (u)≥ 1 is the degree of u.
Note that Φ(X)≤ n·_max.
We let Δ_t=Φ(_t+1)-Φ(_t) be the potential difference in step t.
In addition, let
_t=, and R={ (u,v)∈ E  u∈ and  v∉} be the set of  edges in  with one endpoint being mutant and the other being resident.
Moreover, 
denote
F=∑_u∈ V_ 
as
the total population fitness in .
Given a pair (u,v)∈ R, let p_u→ v be the probability that u reproduces and replaces v.

First we show that (Δ_t)≥ 0, i.e., in expectation, the potential function increases in each step:


    (Δ_t)   =∑_(u,v)∈ R(p_u→ v·(v)/(v) - p_v→ u·(u)/(u))
       =∑_(u,v)∈ R((u)/F1/(u)(v)/(v)-(v)/F1/(v)(u)/(u))
       =∑_(u,v)∈ R(u)((v)-(v))/d(u)d(v)F≥ 0


as (v)≥(v) since  is mutant-biased.
Second, we give a lower bound on the variance of Δ_t when ∅⊂⊂ V, and thus there exists an edge (u,v)∈ R.
First, we have 


    p_v→ u = (v)/F1/(v)≥/n·1/n =/n^2·


while the potential function changes by 
Δ_t≤-(u)/d(u).
Therefore, [Δ_t≤-(u)d(u)]≥/n^2·, and


    (Δ_t)   ≥[Δ_t≤ -(u)/d(u)]·(-(u)/d(u)-(Δ_t))^2
       ≥/n^2·(-/n)^2=·^2/n^4·.


The potential Φ gives rise to a submartingale with upper bound B=n·.
The re-scaled function Φ(Φ-2B)+B^2 satisfies the conditions
of the upper additive drift theorem <cit.> with initial value at most
B^2 and step-wise drift at least (Δ_t). 
We thus arrive at


    T(G,S) ≤B^2/(Δ_t) = n^2·^2/·^2/n^4·≤n^6·^3/^3.




<ref> yields an FPRAS for the fixation probability when mutant and resident fitnesses are polynomially (in n) related.


Given a mutant-biased undirected fitness graph  with /=n^O(1) and a seed set ⊆ V,
the fixation probability _() admits an FPRAS.






§ HARDNESS OF OPTIMIZATION


Here we turn our attention to the seed selection problem, and prove two hardness results.
First, we show that on arbitrary graphs, for any 0<ε<12, it is -hard to distinguish between graphs that achieve maximum fixation probability at most ε and at least 1-ε.
This is in sharp contrast to standard cascade models of influence spread, for which the optimal spread can be efficiently approximated <cit.>.
Then we focus on mutant-biased graphs, and show that achieving the maximum fixation probability remains -hard even in this restricted setting.

Our reduction is from the -hard problem Set Cover <cit.>.
Given an instance (, , ), where  is a universe,  a set of subsets of , and  a size constraint, the task is to decide whether there exist  subsets in  that cover .
Wlog, =⋃_A∈A.
We construct a fitness graph =(G,(,)) where G = (V, E, ) is a bipartite graph with two parts V = V_1 ∪ V_2 with V_1 = and V_2 =, and define the edge relation as 
E = {(u,v) ∈ V_1 × V_2  v ∈ u}∪ (V_2 × V_1)
i.e., there is an edge (u, v) for each element v of  that appears in the set u of , as well as all possible edges from V_2 to V_1.
The weight function is uniform: (u,v)=1/(u) for each (u,v)∈ E.
The resident fitness is (u)=1 for all u∈ V.
The mutant fitness is parametric on two values x≥ 1 and y≤ 1 to be fixed later, with (u)=x if u∈ V_1 and (u)=y if u∈ V_2.
See <ref> for an illustration.

Our construction guarantees upper and lower bounds on the fixation probability depending on whether the seed set forms a set cover of (, ), as stated in the following lemma.

lemmalemhardnessmain
The following assertions hold.


  *  If  is not a set cover then _() ≤ 1-(1n/1n+(n-1)y)^n.

  *  If  is a set cover then


    _() ≥(y/n^2(xn/xn+n)^n/1-(1-y/n^2)(xn/xn+n)^n)^n.




Before we prove <ref>, 
we show 
how <ref> leads to the two hardness results of this section.
theoremthmhardnessgeneral
For any 0<ε< 12, it is -hard to distinguish between 
max__()≤ε and
max__()> 1-ε.


We solve the inequalities of <ref>, 
and obtain that there exist y=1/O(n^3) and x=O(n^10) satisfying them.
As both values are polynomial in n, 
this completes a polynomial reduction from Set Cover to seed selection.




theoremthmhardnessmutantbiased
In the class of mutant-biased fitness graphs, it is -hard to distinguish between max__()≤ 1-1/(n^2n) and max__()> 1-1/(n^2n).


We set y=1 and solve the second inequality of <ref>.
We obtain that it is satisfied by some x=2^O(nlog n).
The fitness graph is mutant-biased as (u)=1 and (u)≥ 1 for all nodes u.
The description of x is polynomially long in n, thus we have a polynomial reduction from Set Cover to seed selection on mutant-biased graphs.


We remark that the class of graphs behind <ref> form a special case of the Positional Moran process <cit.>, by setting as active nodes =V_1 and fitness advantage δ=2^O(nlog n). 
Thus, the -hardness of <ref> extends to the Positional Moran process.

We now turn our attention to the proof of <ref>.
By a small abuse of terminology, we say that a configuration  covers V_2 to denote that the sets in ∩ V_1 cover V_2.
<ref> relies on the following intermediate lemma,
which intuitively states that, starting from a configuration _1 that contains a resident node v∈ V_2 not covered by ,
the process loses all mutants in V_1 with large enough probability.

lemmalemupperresident
From any configuration _1 with V_2∖ (_1 ∪{ u∈_1 (u,v)∈ E })≠∅,
the process reaches a configuration _2 with _2∩ V_1=∅ with probability p≥(1n/1n+(n-1)y)^|V_1|.


We can now prove the upper bound of <ref>.


First, we show that the probability q of reaching configuration _1 such that V_2∖ (_1 ∪{ u∈_1 (u,v)∈ E })≠∅ is at least 1n/1n+(n-1)y.
Then, by using <ref> on _1 we derive that the process reaches a configuration _2 with _2∩ V_1=∅ with probability at least (1n/1n+(n-1)y)^|V_1|=q^|V_1|.
While at configuration _2, the process changes configuration when either a resident in V_1 replaces a mutant in V_2, or vice versa.
Recall that the probability that the first event occurs before the second is at least q.
Repeating the process for all mutants in V_2∖{v} (as v is already a resident in _2) we arrive in a configuration without mutants in V_2 with probability at least q^|V_2|-1.
At this point all mutants have gone extinct, thus
_()≤ 1- q^1+|V_1|+(|V_2|-1)
= 1-(1n/1n+(n-1)y)^n.


The following lemma 
states that, starting from a configuration  that covers V_2, the process makes all nodes in V_2 mutants without losing any mutant in V_1, with certain probability.

lemmalemlowermutant
From any configuration  that covers V_2,
the process reaches a configuration ^* with V_2∪ (∩ V_1)⊆^* with probability p^*≥(xn/xn+n)^|V_2|.


We can now prove the lower bound of <ref>.


We consider 4 configurations; 
any configuration  that covers V_2;
^- with less mutants in V_1 than ;
^* with same mutants in V_1 with  and all nodes in V_2 being mutants;
and
^+ starting from ^* includes at least one more mutant in V_1.
The Markov chain in <ref> captures this process where states S_1, S_2, S_3 and S_4 denote that the process is in configurations ^-, , ^* and ^+, respectively. 
To prove the Lemma, we first bound the transition probabilities of Markov chain in <ref>.
We prove that p^+≥y/n^2(xn/xn+n)^n/1-(1-y/n^2)(xn/xn+n)^n, 
p^*≥(xn/xn+n)^n and q≥y/n^2.
Note that p^+ is lower-bounded by the probability that a random walk starting in S_2 (i.e., ) gets absorbed in S_4 (i.e., ^+). 
Let x_i be the probability that a random walk starting in S_i gets absorbed in S_4. 
We have

x_2 = p^*· x_3 + (1-p^*) · x_1
and x_3 = q · x_4 + (1-q) · x_2,
with boundary conditions x_1 = 0 and x_4 = 1, whence x_2 = q · p^*/1-(1-q) · p^*. 
Since ⊂^+, the set ^+ also covers V_2, thus the reasoning repeats for up to n steps until fixation, resulting in 
_()≥(p^+)^n.








§ MONOTONICITY AND SUBMODULARITY


<ref> rules out polynomial-time algorithms for any non-trivial answer to seed selection.
For mutant-biased graphs, <ref> states that the problem remains -hard to solve exactly, but does not rule out tractable approximations.
Indeed, here we prove that, on mutant-biased graphs, the fixation probability is monotone and submodular, thus seed selection admits a constant-factor approximation.
Our proofs are based on coupling arguments.
Instead of applying these arguments directly to the Heterogeneous Moran process, we introduce a slight variation, the Loopy process, and argue that it is equivalent to the Heterogeneous process in the sense that it preserves the fixation probability.

The Loopy process 
Our Loopy variant is similar to the original process, but we slightly modify the underlying fitness graph =(G,(, )) in each step based on the current configuration . 
Without loss of generality, we let every node u∈ V have a self-loop (u,u)∈ E (we can always assign (u,u) = 0). 
When the original process is at some configuration , different nodes reproduce at different rates. 
When the Loopy process is at configuration , we construct a fitness graph _=(G_, (, )),
where  is the constant function u↦ 1, and G_=(V,E,_) is a graph with the same structure as G, but the weight function is modified by adjusting the self-loop probability of each node, as follows.


    _(u,v)=
    _(u)/_max·(u,v),               if  u≠ v 
    
    1 - _(u)/_max(1-(u,v)),     if  u = v


Intuitively, all 
nodes in _ reproduce at equal rates as
they have the same fitness regardless of their type.
The new weight function _ compensates for this uniformity of reproduction rates: nodes that 
formerly had lower fitness now have stronger self-loops,
which restores the relative propagation rates between neighbors.
That is, the probability distribution [_t+1 |_t+1≠_t],
and thus the fixation probability, is identical in the two processes, as stated in the following lemma.

lemma
lemloopyeq
For any seed set, the Heterogeneous and Loopy Moran processes share the same fixation probability.


Relation to the Two-Graphs process
The Loopy process is a special case of the recent Two-Graphs Moran process <cit.>. 
To obtain the Two-Graphs process,
we define two graphs G_M and G_R for mutants and residents, respectively.
For each edge (u,v), its weight w_M(u,v) in G_M and w_R(u,v) in G_R, is obtained from <ref>, considering that u∈ and u∉, respectively.
In turn, <ref> implies that the hardness of <ref> also hold for seed selection in the Two-Graphs model.

Monotonicity 
The following monotonicity corollary follows from <ref> and the monotonicity of the Two Graphs process <cit.>.

corollary
lemmonotonicity
For any fitness graph =(G,(, )) and any two seed sets ⊆', we have _()≤_(').



Submodularity We now turn our attention to the submodularity of the fixation probability in the Heterogeneous Moran process.
Although the function is not submodular in general, we prove that it becomes submodular on mutant-biased fitness graphs.
In particular, we show that for any two seed sets , ⊆ V, the following submodularity condition holds:


    _() + _() ≥_(∪) + _(∩)


Our proof is via a four-way coupling of the corresponding processes starting in one of the seed sets of <ref>.

lemmalemsubmodularity
For any mutant-biased fitness graph =(G,(, )), the fixation probability _() is submodular.


Let _1 = (^1_t)_t ≥ 0, _2 = (^2_t)_t ≥ 0, _3 = (^3_t)_t ≥ 0, and _4 = (^4_t)_t ≥ 0, be four Loopy processes with seed sets , , ∪ and ∩, respectively. 
To prove submodularity, we employ two tricks for _3. 
First, along its configurations ^3_t, we also keep track of the set of mutants _t (resp., _t) that are copies of some initial node in  (resp., ). 
Whenever a node v receives the mutant trait from a neighbor u, we place v in _t+1 (resp., _t + 1) following the membership of u in _t (resp., _t). 
Initially, _0 = and _0 =. 
Second, with probability 1, every run of _3 that results in fixation, eventually (i.e., if we let the process run on) leads to the fixation of  or  
(possibly both, assuming ∩≠∅); 
that is, every node is a copy of some node in  or . 
We thus compute the fixation probability with seed ∪ by summing over runs in which  or  fixates.

To prove submodularity, we consider this refined view of the process and establish a four-way coupling between _1, _2, _3 and _4 that guarantees the following invariants: (i) ^1_t∪^2_t⊆^3_t, (ii) ^4_t⊆^1_t∩^2_t, (iii) _t⊆^1_t, and (iv) _t⊆^2_t. Now, consider any execution in which _3 fixates. Since  or  eventually fixates in _3, due to invariants (iii) and (iv), at least one of _1, _2 fixates as well. Moreover, if _4 also fixates, due to invariant (ii),
both _1 and _2 fixate. Thus the invariants guarantee submodularity.

The invariants hold at t=0. 
Now, consider some arbitrary time t with the four processes at configurations ^j_t = ^j, for j∈{1,2,3,4}, _t =, and _t =. 
To obtain ^j_t + 1, we sample the same node u for reproduction with probability 1n in all processes. 
From invariants (i) and (ii), and since (u)≥(u), we derive that __3(u,u) ≤__j(u,u) ≤__4(u,u) for j ∈{1, 2}, as residents have a larger self-loop weight.
In _3, we choose a neighbor v of u with probability __3(u, v) and propagate the trait of u to v. 
In _1, _2 and _4, if u = v, we perform the same update; otherwise, if u has the same type as in _3, we also perform the same update. 
From the invariants, if u is resident in _3 then the same holds in all other processes, while if u is a mutant in _3 then the same holds in at least one of _1, _2 (depending on whether u ∈ and u ∈), and if that holds for _1 and _2, then it holds for _4. 
However, if u is resident in _j for some j∈{1, 2} but mutant in _3, i.e., u ∈^3 ∖^j, then, due to invariant (ii), u is also resident in _4, i.e., u ∈^3 ∖^4; then, in _j and _4, u propagates to itself with probability __j(u,u)-__3(u,u)≥ 0, 
and to v with the remaining probability 1 - (__j(u,u) - __3(u,u)). 
It follows that all three invariants are maintained.


Following <cit.>, monotonicity and submodularit lead to the following approximation guarantee.


Given a mutant-biased 
fitness graph  and budget , let ^* be an optimal seed set and _gr the solution of the Greedy algorithm. 
We have _ (_gr) ≥( 1 - 1e) _ (^*).


The Greedy algorithm builds the seed set iteratively by choosing the node that yields the maximum fixation probability gain.
Finally, note that due to symmetry, on resident-biased fitness graphs  ((u)≤(u) for all u), 
_() is supermodular, thus Greedy offers no approximation guarantees.









§ EXPERIMENTAL ANALYSIS


Here, we present our experimental evaluation of the Greedy algorithm and other network heuristics, 
varying the seed size k and the maximum mutant fitness m_max.




Datasets We use four real-world networks from  Netzschleuder, SNAP and Network Repository (<ref>).



  * Facebook: A Facebook ego network in which nodes represent profiles and edges indicate friendship.

  * Colocation: A proximity network of students and teachers of a French school.
Edge weights count
the frequency of contact between individuals
during a two-day period.

  * Mammalia: An animal-contact network based on movements of voles (Microtus agrestis).
Each edge weight counts the common traps the two voles were caught in.

  * Polblogs:
A network of hyperlinks among a large set of U.S. political weblogs from before the 2004 election. 

Our experiments are not meant to be exhaustive, but rather indicative of the performance of the greedy algorithm and common network-optimization heuristics on a few diverse networks.
We set the resident fitness to 1,
while the mutant fitness of each node u is determined by sampling a uniform distribution (u)∼𝒰(1,).
This results in mutant-biased graphs, for which <ref> guarantees that the fixation probability admits a Monte Carlo approximation.

Greedy and Baselines
We evaluate the performance of the standard Greedy algorithm behind <ref> <cit.> 
against four common baseline algorithms from related literature on seed selection under diffusion processes <cit.>.


  * Random: select uniformly at random.

  * Degree: select by smallest degree.

  * Closeness: select by smallest closeness centrality.

  * PageRank: select by smallest PageRank score.



The Random selection strategy is a standard baseline to measure the intricacy of the problem.
Degree is the only existing algorithm for seed selection in the Moran model, and is optimal for undirected and unweighted networks under the neutral setting
(but underperforms when >1).
On the other hand, Closeness and PageRank take into account the structure of the graph and its connectivity.
For these two centrality heuristics we also tried selecting the top-k-nodes by largest value, which resulted in worse performance.
All Monte Carlo simulations were run over 5000 iterations.

Performance vs.  <ref> shows performance as the size constraint k increases for a fixed mutant fitness distribution. 
In agreement with <ref>, the performance of all algorithms rises as k grows, while Greedy has diminishing returns. 
Notably, Greedy outperforms all heuristics especially for small size constraints, while  Pagerank forms high quality solutions for the undirected and unweighted graph Facebook. 
On the other hand, seed selection becomes more challenging for directed (Mammalia, Polblogs) and edge-weighted graphs (Colocation, Mammalia), in which only Greedy uncovers high-quality seed sets.


Performance vs.  <ref> shows performance as the mutant fitness interval [1,] increases, for fixed size .
Random selection performs poorly, showing that the problem is not trivial,
while the other two heuristics have mixed performance.
On the other hand, Greedy achieves
a steady, high-quality performance in all datasets and problem parameters.





§ CONCLUSION

We studied a natural optimization problem pertaining to network diffusion by the Heterogeneous Moran process, namely selecting a set of seed nodes that maximize the effect of the invasion.
To our knowledge, this is the first paper to study this standard optimization problem on Moran models.
We showed that the problem is strongly inapproximable in general, but becomes 
approximable 
on mutant-biased graphs, although 
the exact solution remains -hard.
Several interesting questions remain open for future work, such as,
is seed selection hard in the Standard model; 
and are there tighter approximations for mutant-biased graphs?




§ ACKNOWLEDGMENTS

A.P. was partially supported by a research grant (VIL42117) from VILLUM FONDEN.
J.T. was supported by Charles Univ. projects UNCE 24/SCI/008 and PRIMUS 24/SCI/012.
named







§ APPENDIX


First we prove <ref> and <ref>, that are used in the analysis of <ref>.

*

Let q=1n/1n+(n-1)y.
While at _1, there exists at least one resident node v∈ V_2 that is not covered by any mutant in V_1.
Let '_1 be the first configuration that the process reaches in which the number of mutants in V_1 has changed from _1
(i.e., V_1 either has one more, or one less mutant in '_1 compared to _1).
Let F=∑_u∈ V__1(u) be the total population fitness at _1.
The probability that v replaces a mutant in V_1 in a single step is p_1≥1/F1/n.
On the other hand, the probability that any mutant in V_2 replaces a resident in V_1 in a single step is p_2≤|V_2| y/F≤(n-1) y/F.
Thus, the probability that V_1 has lost a mutant in '_1 is at least p_1/p_1+p_2≥1n/1n+(n-1)y=q.

Now, observe that '_1 satisfies the conditions of _1, i.e., V_2∖ ('_1 ∪{ u∈'_1 (u,v)∈ E })≠∅.
Thus we can repeat the above process until all nodes in V_1 have become residents, leading to the desired configuration _2 with probability


    p≥ q^|V_1|=(1n/1n+(n-1)y)^|V_1|





*

Let q=x/n/x/n+n, and let F=∑_u∈ V_(u) be the total population fitness in .
As long as there are residents in V_2, the probability that a mutant node in V_1 replaces a resident in V_2 is at least p_1≥x/F1/n.
On the other hand, the probability that a resident in V_2 (resp. V_1)
replaces a mutant in V_1 (resp. V_2) is at most p_2≤|V_2|/F (resp. p_3≤|V_1|/F).
The probability that the first event occurs before the second and third one is thus at least p_1/p_1+p_2+p_3≥x/n/x/n+|V_2|+|V_1|=x/n/x/n+n=q.
Observe that any other event that changes the configuration (i.e., mutants in V_2 replacing residents in V_1) 
results in a configuration that also covers V_2, 
thus we can repeat the above argument until all nodes in V_2 become mutants, which occurs with probability at least


    p^*≥ q^|V_2|≥(xn/xn+n)^|V_2|






 We continue with proving the upper and lower bounds of  <ref>.
*

We prove the two assertions separately.

    
  * Item 1:
    First, we argue that with probability at least q=1n/1n+(n-1)y, the process reaches a configuration _1 such that  V_2∖ (_1 ∪{ u∈_1 (u,v)∈ E })≠∅.
    Indeed, since  does not form a set cover, there exists a node v∈ V_2 that has no mutant incoming neighbor in V_1.
    If v∉, we are done.
    Otherwise, the probability that, in a single step, v is replaced by any resident neighbor in V_1 is at least p_1≥1/F1/n, where F=∑_u∈ V_(u) is the total population fitness at .
    On the other hand, any resident in V_1 is replaced by mutants in V_2 with probability p_2≤|V_2|y/F≤(n-1)y/F.
    Thus, the probability that the process reaches a desired configuration _1 is at least p_1/p_1+p_2≥1/n/1/n+(n-1)y=q.
    
    Second, <ref> applies on _1 to show that the process reaches a configuration _2 with _2∩ V_1=∅ with probability at least (1n/1n+(n-1)y)^|V_1|=q^|V_1|.
    
    Third, while at configuration _2, the process changes configuration when either a resident in V_1 replaces a mutant in V_2, or vice versa.
    We have already argued in the first step that the probability that the first event occurs before the second is at least q.
    Now, we repeat this process until all mutants in V_2 have become residents, which occurs with probability at least q^|V_2|-1, as v∈ V_2 is already a resident in _2.
    At this point the mutants have gone extinct, thus
    _()≤ 1- q^1+|V_1|+(|V_2|-1)
    = 1-(1n/1n+(n-1)y)^n.

    
  * Item 2: We first prove the following statement. 
    Consider the process at any configuration  that covers V_2, and let p^+ be the probability that it reaches a subsequent configuration ^+ with at least one more mutant in V_1. 
    We will show that p^+≥y/n^2(xn/xn+n)^n/1-(1-y/n^2)(xn/xn+n)^n.
    Given any such configuration , let p^* be the probability that the process reaches a subsequent configuration ^* with V_2 ∪ (∩ V_1)⊆^*. 
    By <ref>, we have p^*≥(xn/xn+n)^|V_2|≥(xn/xn+n)^n.
    While at ^*, 
    the process only progresses 
    when a mutant in V_2 replaces a resident in V_1, or a resident in V_1 replaces a mutant in V_2. 
    As long as there are residents in V_1, the probability q that a mutant from V_2 replaces a resident in V_1 before any such resident reproduces satisfies
    q≥y/F1/n/y/F1/n+|V_1|/F≥y/F1/n/y/F1/n+n-1/F=y/y+(n-1)n≥y/n^2, where the last inequality holds as y≤ 1≤ n. 
    If this happens, we reach the desired configuration ^+.
    Otherwise, 
    a resident in V_1 replaces a mutant in V_2, the resulting configuration still covers V_2, and the argument repeats. 
    The Markov chain in <ref> captures this process. 
    The states S_2, S_3 and S_4 denote that the process is in configurations , ^* and ^+, respectively. 
    Hence, p^+ is lower-bounded by the probability that a random walk starting in S_2 (i.e., ) gets absorbed in S_4 (i.e., ^+). 
    Let x_i be the probability that a random walk starting in S_i gets absorbed in S_4. 
    We have
    x_2 = p^*· x_3 + (1-p^*) · x_1
    and x_3 = q · x_4 + (1-q) · x_2,
    with boundary conditions x_1 = 0 and x_4 = 1, whence x_2 = q · p^*/1-(1-q) · p^*. 
    Since ⊂^+, set ^+ also covers V_2, thus the reasoning applies for up to n steps until fixation, resulting in 
    _()≥(p^+)^n≥(y/n^2(xn/xn+n)^n/1-(1-y/n^2)(xn/xn+n)^n)^n.




We continue with the proof of <ref> and <ref>.
To make our analysis easier, we first prove two simple lemmas.

lemmaloglowerbound
For every ζ>0, we have ln(1+1/ζ)≥1/ζ+1.




    (1+1/ζ)^ζ+1≥ e
    ⇔   ln(1+1/ζ)^ζ+1≥ln e
    ⇔   ln(1+1/ζ)≥1/ζ+1






lemmalemprobbound
Let p∈(0,1), and β≤ln(1p)/n.
We have


    ( 1/1+β)^n≥ p




Let γ=1+1β, thus β=1(γ-1).
Then 


    γ = 1 + 1/β≥ 1+ n/ln(1/p)


Moreover


    1/1+β = 1/1+1/γ-1=1/γ/γ-1=γ-1/γ=1-1/γ


Let k=ln(1p)^-1, thus p=e^-1k and γ≥ nk+1.
We have


    (1-1/nk + 1)^nk≥ e^-1
    ⇒   (1-1/nk + 1)^n ≥ e^-1/k
    ⇒   ( 1-1/γ)^n≥ p 
    ⇒   ( 1/1+β)^n ≥ p


as desired.



*

The proof  is by algebraic manipulation on the inequalities of <ref>.
In particular, we argue that there exist x and y that have polynomially-long description (in n) for which the inequalities stated in the theorem hold.
In turn, this completes a polynomial reduction from Set Cover to distinguishing between max__()≤ε and
max__()> 1-ε in the Heterogeneous Moran process.

First, assume that  does not form a set cover, thus <ref> of <ref> applies.
We solve the corresponding inequality to arrive at a suitable value for y.
In particular, for _()≤ε, it suffices to find a y small enough such that


    1-(1/n/1/n+(n-1)y)^n≤ε
    ⇔    1-ε≤(1/1+n(n-1)y)^n


Using <ref> for p=1-ε, we set


    n(n-1)y =β≤ln(1/p)/n
    ⇒    y≤1/O(n^3)


as ε is fixed.
Hence y suffices to be polynomially small in n for the inequality of <ref> of  <ref> to hold.


On the other hand, if  forms a set cover, <ref> of <ref> applies.
We solve the corresponding inequality to arrive at a suitable value for x.
In particular, for _() >1-ε, it suffices to find an x large enough such that


    (y/n^2(x/n/x/n+n)^n/1-(1-y/n^2)(x/n/x/n+n)^n)^n >1-ε


Substituting with α=(x/x+n^2)^n, we have


    (y/n^2α/1-(1-y/n^2)α)^n >1-ε
    ⇒    ( y/n^2α/1-α +y/n^2α)^n >1-ε
    ⇒    ( 1/1+ n^2(1-α)/yα)^n > 1-ε


Using <ref> for p=1-ε, we set


    n^2(1-α)/yα=β≤ln(1/p)/n
    ⇒    α≥n^3/n^3+yc


where c=ln(1p) is a constant.
We thus have


    (x/x+n^2)^n≥n^3/n^3+yc
    ⇒    (1/1+n^2/x)^n ≥n^3/n^3+yc


Using <ref> for p=n^3/n^3+yc, we set


    n^2/x = β≤ln(1/q)/n
    ⇒    x≥n^3/ln(1/q) = n^3/ln(1+1/n^3/cy)


Using <ref> for ζ=n^3/cy=O(n^6) as c=O(1) and y=1/O(n^3), we have


    x≥n^3/1/n^7=O(n^10)


Thus x suffices to be polynomially large in n for the inequality of <ref> of <ref> to hold.


*

The proof  is by algebraic manipulation on the inequalities of <ref> for the specific case where y=1 and x≥ y.
Observe that this makes the corresponding graph mutant-biased.
In particular, we argue that there exists an x with a polynomially-long description (in n) for which the inequalities stated in the theorem hold.
In turn, this completes a polynomial reduction from Set Cover to distinguishing between max__()≤ 1-1/(n^2n) and max__()> 1-1/(n^2n) in the Heterogeneous Moran process.

First, assume that  does not form a set cover.
<ref> of <ref> for y=1 gives


    _() ≤ 1-(1/1+(n-1)n)^n≤ 1-1/n^2n



On the other hand, if  forms a set cover, <ref> of <ref> applies.
For y=1, we solve the corresponding inequality to arrive at a suitable value for x.
In particular, for _() >1-1/n^2n, it suffices to find an x large enough such that


    _() ≥(1/n^2(x/n/x/n+n)^n/1-(1-1/n^2)(x/n/x/n+n)^n)^n>1-1/n^2n


Substituting with a=(x/x+n^2)^n, we have


    (1/n^2α/1-(1-1/n^2)α)^n >1-1/n^2n
    ⇒    ( 1/n^2α/1-α +1/n^2α)^n >1-1/n^2n
    ⇒    ( 1/1+ n^2(1-α)/α)^n > 1-1/n^2n



Using <ref> for p=1-1/n^2n, we set


    n^2(1-α)/α=β≤ln(1/p)/n
    ⇒    α≥n^3/n^3+c


where c=ln(1p).
We thus have


    (x/x+n^2)^n≥n^3/n^3+c
    ⇒    (1/1+n^2/x)^n ≥n^3/n^3+c



Using <ref> for p=n^3/n^3+c, we set


    n^2/x = β≤ln(1/q)/n
    ⇒    x≥n^3/ln(1/q) = n^3/ln(1+1/n^3/c)



Using <ref> for ζ=n^3/c, we have



    x≥n^3/1/n^4/c=n^7/c



Finally, recall that c=ln(1p)=ln(1+1/n^2n-1).
Using <ref> again with β=n^2n-1, we conclude that



    x≥n^7/ln(1+1/n^2n-1)≥n^7/1/n^2n=2^O(nlog n)


Observe that x has polynomially-long (in n) description, thus our reduction from Set Cover is in polynomial time.


*

Consider the Heterogeneous and Positional Moran process, and assume that they are in the same configuration ⊆ V.
Consider any edge (u,v)∈ E with u≠ v.
Let p_u→ v and p'_u→ v be the probabilities of u transferring its trait to v under the Heterogeneous and Loopy Moran processes, respectively, when in configuration . 
By the definition of the models, we have


    p_u→ v = _(u)/w(u,v),     p'_u→ v = 1/nw_(u,v)=_(u)·(u,v)/n·_max.


Moreover, let R=E∖{(x,x) x∈ V} be the set of edges without the self loops in G.
Note that, form , each process can progress to a distinct configuration '≠ only if a node u transfers its trait along an edge (u,v)∈ R.
Let p_1 and p'_1 be the probability that this occurs in the Heterogeneous and Loopy process, repsectively, and we have


    p_1=∑_(x,y)∈ R_(x)/w(x,y),          p_1'=   ∑_(x,y)∈ R1/nw_(x,y)
    =   ∑_(x,y)∈ R_(x)·(x,y)/n·_max.



Finally, observe that


    p_u→ v/p_1=_(u)/w(u,v)/∑_(l,r)∈ R_(l)/w(l,r)=_(u)w(u,v)/∑_(l,r)∈ R_(l)w(l,r)
    p'_u→ v/p'_1=_(u)·(u,v)/n·_max/∑_(l,r)∈ R_(l)·(l,r)/n·_max=_(u)w(u,v)/∑_(l,r)∈ R_(l)w(l,r).



Thus, the probability distribution [_t+1 = |_t+1≠_t] is the same in the two processes, yielding the same fixation probability starting from the same seed set, as desired.

