Article abstract
Article abstract
Adaptive evolution of dispersal strategies is one mechanism by which species can respond to rapid environmental changes. However, under rapid anthropogenic fragmentation, the evolution of dispersal may be limited, and species may be unable to adequately adapt to fragmented landscapes. Here, we develop a spatially explicit model to investigate the evolution of dispersal kernels under various combinations of fragmentation dynamics and initial conditions. We also study the consequences of modelling an evolutionary process in which dispersal phenotypes continuously and gradually shift in phenotype space in a manner corresponding to a polygenic underlying genetic architecture. With rapid fragmentation rates, we observed the emergence of long-term transient states in which dispersal strategies are not well suited to fragmented landscapes. We also show that the extent and length of these transient states depend on the pre-fragmentation dispersal strategy of the species, as well as on the rate of the fragmentation process leading to the fragmented landscape. In an increasingly fragmented world, understanding the ability of populations to adapt, and the effects that rapid fragmentation has on the evolution of dispersal, is critical for an informed assessment of species viability in the Anthropocene.
Introduction
Introduction
This notebook includes code for running the simulations detailed in the article. The “RunOneSim” function summarizes all the different parts of the simulation, and returns the results of a single evolutionary iteration of the process (the evolution of dispersal in a population undergoing fragmentation). To obtain the results described in the manuscript one should run multiple repeats of the process and summarize the results across the simulation replicates. in addition, at the end of the notebook are three visualization plots of the results: (1) The mean dispersal kernel of the population at the end of the simulation, (2) The mean proportion of non-dispersers (the fraction of the distance class 0 in the dispersal kernels) over time, and (3) the seed mortality (proportion of dispersed seeds reaching uninhabitable cells) over time.
Code for running a single iteration of the process of evolution of dispersal in a landscape undergoing dynamic fragmentation
Code for running a single iteration of the process of evolution of dispersal in a landscape undergoing dynamic fragmentation
The code below describes the different elements that go into the “RunOneSim” function. This includes: (1) generation of a random autocorrelated landscape, (2) initiation of the population with one of 4 possible initial kernels, (3) simulation of the death-and-replacement process through dispersal, (4) mutation of the dispersal kernel, and (5) a dynamic landscape fragmentation process.
$HistoryLength=0;
Initiate age at 1-125 at each habitable cell
IntiateAgesMatrix[landmat_]:=Module[{mat},mat=landmat/.1RandomInteger[{1,125}];mat]
CreateRandomLandscape[n_,Phabitat_]:=Module[{a1,m1,m2,m3},a1=Sort[RandomSample[Flatten[Table[{i,j},{i,n},{j,n}],1],Round[Phabitat]]];m1=SparseArray[a1Table[1,{Length[a1]}]];m2=Normal[m1];m2]
2
n
Nextevent[ladscapematAges_,habitatpos_]:=Module[{ages,t1,t2,Minage},ages=Table[ladscapematAges〚k〚1〛,k〚2〛〛,{k,habitatpos}];Minage=Min[ages];t1=Position[ladscapematAges,x_/;xMinage];t2=t1⋂habitatpos;RandomChoice[t2]]
Initial kernel - Uniform
InitiateDispK0[ladscapematPlain_,s_]:=Module{w1},w1=ladscapematPlain/.1;w1
Table[1,{s}]
s
s is the number of bins in the dispersalK histogram. We now initiate everything with uniform kernels. Kernels must sum to 1
Initial kernel - 1 in bin 0, 0 otherwise (no dispersal)
InitiateDispK1[ladscapematPlain_,s_]:=Module[{w1},w1=ladscapematPlain/.1Table[If[i1,1,0],{i,s}];w1]
Initial kernel - 0.5 in bin 0 and s-1, 0 otherwise (bimodal dispersal)
InitiateDispK2[ladscapematPlain_,s_]:=Module[{w1},w1=ladscapematPlain/.1Table[If[i1||is,0.5,0],{i,s}];w1]
Initial kernel - 1 in last bin, 0 otherwise (long-distance dispersal)
InitiateDispK3[ladscapematPlain_,s_]:=Module[{w1},w1=ladscapematPlain/.1Table[If[is,1,0],{i,s}];w1]
TorusEucludianDistance[a_,b_,n_]:=Module{xdist,ydist},xdist=If[a〚1〛≤b〚1〛,Min[{b〚1〛-a〚1〛,a〚1〛+n-b〚1〛}],Min[{a〚1〛-b〚1〛,b〚1〛+n-a〚1〛}]];ydist=If[a〚2〛≤b〚2〛,Min[{b〚2〛-a〚2〛,a〚2〛+n-b〚2〛}],Min[{a〚2〛-b〚2〛,b〚2〛+n-a〚2〛}]];+
2
xdist
2
ydist
Select a new tree position that will replace a dying tree (in Treeposition)
CollectSeeds[Treeposition_,dispKmat_,habitatposList_,s_]:=Module{potentialdoners},potentialdoners=Select[habitatposList,Floor[TorusEucludianDistance[Treeposition,#1,Length[dispKmat]]]≤s-1&];seedpile=TableIfdispKmatt1,t2===0,0,IfTreeposition===t,1,dispKmatt1,t2,Floor[TorusEucludianDistance[Treeposition,t,Length[dispKmat]]]+1,{t,potentialdoners};If[Total[seedpile]===0||Total[seedpile]===0.,seedpile=Table[1,Length[seedpile]]];RandomChoice[seedpilepotentialdoners]
1
2πTorusEucludianDistance[Treeposition,t,Length[dispKmat]]
MutateDispK[dispK_,mu_,magnitude_]:=Module{t1},t1=TableIfRandomReal[]<,Max[RandomReal[{v-magnitude,v+magnitude}],0],v,{v,dispK};IfTotal[t1]===0||Total[t1]===0.,Table,{Length[dispK]},(*Normalizesdispersalkerneltointegrateto1*)
mu
Length[dispK]
1
Length[dispK]
t1
Total[t1]
Here each bin in the dispersal kernel is tested if it mutates, then how much, and at the end the kernel s normalized.
mu is the probability of a mutation occurring in one bin in one generation (more precisely, since it can be larger than one - the number of bins expected to change each generation); Magnitude is the maximal size of the mutation
mu is the probability of a mutation occurring in one bin in one generation (more precisely, since it can be larger than one - the number of bins expected to change each generation); Magnitude is the maximal size of the mutation
DeathStep[Treeposition_,dispKmat_,habitatposList_,landdata_,s_,AgeMatrix_,mu_,magnitude_,time_]:=Module[{newtime,newseed,newAgeMatrix,newdispK,newdispKmat},(*Conductthereplacemntofa=nindividualbyanotherindividual*)newtime=time+AgeMatrix〚Treeposition〚1〛,Treeposition〚2〛〛;newseed=CollectSeeds[Treeposition,dispKmat,habitatposList,s];newAgeMatrix=(AgeMatrix-AgeMatrix〚Treeposition〚1〛,Treeposition〚2〛〛)landdata;newAgeMatrix=ReplacePart[newAgeMatrix,TreepositionRandomReal[NormalDistribution[100,5]]];newdispK=MutateDispK[dispKmat〚newseed〚1〛,newseed〚2〛〛,mu,magnitude];newdispKmat=ReplacePart[dispKmat,TreepositionIf[newseed===0,0,newdispK]];{newAgeMatrix,newdispKmat,newtime}]
Generate autocorrelated landscape
landscapeShape[Hu_,size_]:=Module{r,a,phi,inq,rinq},q=ConstantArray[0,{size,size}];DoDog={If[i0,0,size-i],If[j0,0,size-j]};r=RandomReal[NormalDistribution[]];a=Ifi0&&j0,0,r(-Hu-1);phi=RandomReal[{0,2π}];o1=aCos[phi]+aSin[phi];o2=aCos[phi]-aSin[phi];If[i0&&j0,o1=Im[o1]];Ifi0||i&&j0||j,o1=Re[o1];If[q〚i+1,j+1〛0,q=ReplacePart[q,{i+1,j+1}o1];q=ReplacePart[q,g+1o2]];,{j,0,size-1},{i,0,size-1};inq=InverseFourier[q];rinq=Re[inq]
1
2
+
2
Min[{i,g〚1〛}]
2
Min[{j,g〚2〛}]
size
2
size
2
Generate random landscape without autocorrelation
landscaperandom[size_]:=Module[{},Table[Table[Random[],size],size]]
Fragmentation rate is measured in proportion of landscape lost per time-step; final frag is the final fragmentation rate, i.e. the p wherein fragmentation stops*)Module[{p,thr,m1},(*Creates landscape of size*size with Hurst index Hu and with p proportion of occupied cells
landscapeChanging[landscapeShape_,time_,fragrate_,finalfrag_,size_]:=Module[{p,thr,m1},(*Createslandscapeofsize*sizewithHurstindexHuandwithpproportionofoccupiedcells*)p=Max[{1-timefragrate,finalfrag}];If[0<p<1,thr=Quantile[Flatten[landscapeShape],1-p];m1=Map[If[#1≤thr,0,1]&,landscapeShape,{2}]];If[p0,m1=ConstantArray[0,{size,size}]];If[p1,m1=ConstantArray[1,{size,size}]];m1]
Number of cells in distance i from a focal cell
gridsatdist[r_]:=Module[{ta2,ta1},ta1=Flatten[Table[{i,j},{i,-r,r},{j,-r,r}],1];ta2=Length[Select[ta1,Floor[EuclideanDistance[{0,0},#1]]r&]];ta2]
habitaatdist[Treeposition_,r_,habitatposList_,dispKmat_,Ngirdastdistlist_]:=Module[{ta1,ta2,ta0},ta0=Length[Select[habitatposList,Floor[TorusEucludianDistance[Treeposition,#,Length[dispKmat]]]==r&]];N[ta0/Ngirdastdistlist[[r+1]]]]
lostseedspos[Treeposition_,habitatposList_,dispKmat_,s_,Ngirdastdistlist_]:=Module[{t1},t1=Table[dispKmat〚Treeposition〚1〛,Treeposition〚2〛,r+1〛habitaatdist[Treeposition,r,habitatposList,dispKmat,Ngirdastdistlist],{r,0,s-1}];1-Total[t1]]
lostseeds[habitatposList_,dispKmat_,s_,Ngirdastdistlist_]:=Module[{},Mean[Table[lostseedspos[i,habitatposList,dispKmat,s,Ngirdastdistlist],{i,habitatposList}]]]
RunOneSim[initialkernel_,timelimit_,s_,mu_,mag_,Hu_,size_,fragrate_,finalfrag_]:=Module[{sims,meanDispKernelList,Ngirdastdistlist,time,landShape,landdata,agemat,habitatpos,DispKmat,flag1,flag2,step,mortalitylist,meanDispKmat,kernelat0,timelist1,timelist2,timelist3,ne,ds,OccupiedDispKernel,res},sims={};Ngirdastdistlist=Table[gridsatdist[r],{r,0,s-1}];time=0;landShape=If[Hu=="Random",landscaperandom[size],landscapeShape[Hu,size]];landdata=landscapeChanging[landShape,0,fragrate,finalfrag,size];(*Intiallandscape-uniformsincetime=0*)agemat=IntiateAgesMatrix[landdata];habitatpos=Sort[Position[landdata,1]];DispKmat=Which[initialkernel==0,InitiateDispK0[landdata,s],initialkernel==1,InitiateDispK1[landdata,s],initialkernel==2,InitiateDispK2[landdata,s],initialkernel==3,InitiateDispK3[landdata,s]];flag1=0;flag2=0;step=1;mortalitylist={lostseeds[habitatpos,DispKmat,s,Ngirdastdistlist]};meanDispKmat=Mean[DeleteCases[Flatten[landdata*DispKmat,1],0]];meanDispKernelList={meanDispKmat};kernelat0={meanDispKmat[[1]]};timelist1={0};timelist2={0};timelist3={0};While[time<timelimit,ne=Nextevent[agemat,habitatpos];ds=DeathStep[ne,DispKmat,habitatpos,landdata,s,agemat,mu,mag,time];agemat=ds[[1]];DispKmat=ds[[2]];time=ds[[3]];If[time>=flag1,OccupiedDispKernel=Select[Flatten[DispKmat,1],Total[#]==1.&];meanDispKmat=Mean[OccupiedDispKernel];AppendTo[meanDispKernelList,meanDispKmat];AppendTo[kernelat0,meanDispKmat[[1]]];AppendTo[timelist2,time];flag1=flag1+100;];If[time>=flag2,AppendTo[mortalitylist,lostseeds[habitatpos,DispKmat,s,Ngirdastdistlist]];AppendTo[timelist3,time];flag2=flag2+1000];landdata=landscapeChanging[landShape,time,fragrate,finalfrag,size];habitatpos=Sort[Position[landdata,1]];agemat=agemat*landdata;DispKmat=Table[landdata[[i,j]]*DispKmat[[i,j]],{i,Length[landdata]},{j,Length[landdata]}]];meanDispKmat=Mean[DeleteCases[Flatten[landdata*ds[[2]],1],0]];res={meanDispKernelList,kernelat0,mortalitylist};res]
The cell below runs a single evolutionary trajectory based on the parameters in the green cell above
Final evolved dispersal kernel (from a single trajectory)
Final evolved dispersal kernel (from a single trajectory)
Proportion of non-dispersers over time (a single evolutionary trajectory)
Proportion of non-dispersers over time (a single evolutionary trajectory)
Dispersal-related seed mortality over time (a single evolutionary trajectory)
Dispersal-related seed mortality over time (a single evolutionary trajectory)