Logistic Map and its Bifurcation Diagram

Introduction

The logistic map is a very simple quadratic map, which also is capable of highly chaotic behavior. It’s equation is like the following:

(Subscript[x, n+1])=r (Subscript[x, n])(1-(Subscript[x, n]))

The constant r is also known as biotic potential, from the fact that this equation is also used as the idealistic population growth model in the field of theoretical ecology. Being idealistic, the model only considers intraspecific interaction based on population density and assumes generations don’t overlap.

Logistic Map is Chaotic

The logistic map is chaotic. Its behaviour significantly varies depending on the value of r. For 0 < r < 1, the population will die eventually, but, for 1 < r ≤ 3, the population will converge to ((r-1)/r). Beyond 3 ≤ r < 4, it will not converge mostly, and further beyond r ≥ 4, it will diverge.

Let’s see what exactly happens:

Module[{x,n},
LogisticMapSequence[r_,len_,x0_]:=RecurrenceTable[{
x[n+1]==r x[n] (1-x[n]),
x[0]==x0
},x,{n,1,len}]
]

ParallelTable[
ListPlot[LogisticMapSequence[r,40,0.88],
Joined->True,
PlotRange->{{1,40},{0,1}},
PlotLabel->HoldForm[r]==r
],
{r,{0.6,1.3,2.95,3.3,3.5,4}}
]//Partition[#,3]&//GraphicsGrid[#,ImageSize->Full]&

For r = 0.6, the population died quickly. For r = 1.3, the population converges quickly, but, for r = 2.95, converges slowly. For r = 3.3, you can observe the oscillation between 2 values, and between 4 values for r = 3.5. Finally, for r = 4, there’s just no pattern.

Since the population diverges for r ≥ 4, a different approach should be used to clearly see what happens. The below plots (n, log|(Subscript[x, n])|):

Log[Abs[#]*10]&/@LogisticMapSequence[4.1,33,0.4]//
ListPlot[#,Joined->True,PlotRange->All]&

The Bifurcation Diagram of Logistics Map

While it’s certainly not impossible to poke into the details of logistics map with the above approach, we certainly need a better view on the subject. This is where “Bifurcation Diagram” kicks in. From Wikipedia:

a bifurcation diagram shows the values visited or approached asymptotically 
(fixed points, periodic orbits, or chaotic attractors) of a system

In simpler words, given a system F, a bifurcation diagram plots all possible values of F(x; σ) for bifurcation parameter x and initial condition σ. That is, to plot the bifurcation diagram of the logistic map, we have to

<<CompiledFunctionTools`
<<Developer`

imageWidth=2560.*4;
imageHeight=1080.*4;
rStart=3.5;
rEnd=4.0;
rLength=rEnd-rStart;
rStep=(rLength/imageHeight)/4.;
seriesSkip=1000;
seriesDepth=Round[.3×imageWidth];

f=Compile[{{x,_Real},{r,_Real}},
r x(1.-x)
];

The following function returns the population values of n-th generations, where skip < n ≤ skip + depth.

runSeries=Compile[{{x0,_Real},{r,_Real},{skip,_Integer},{depth,_Integer}},
Module[{x},
x=Nest[f[#,r]&,x0,skip];
Table[{r,x=f[x,r]},{i,1,depth}]
],
CompilationOptions->{"InlineExternalDefinitions"->True}
];

pointsPerR=ParallelTable[
runSeries[0.5,r,seriesSkip,seriesDepth],
{r,rStart,rEnd,rStep}];
points=Join@@pointsPerR;

plotPoints=Compile[{{pts,_Real,2},{imgW,_Integer},{imgH,_Integer}},
Module[{image,ptx,pty},
With[{rStart=rStart,rEnd=rEnd},
image=Table[0,{imgH+2},{imgW+2}];
Do[
ptx=Round[(pt[[1]]-rStart)/(rEnd-rStart)*imgW]+1;
pty=Round[pt[[2]]*imgH]+1;
image[[pty,ptx]]+=1;,
{pt,pts}
];
image
]],
CompilationOptions->{"InlineExternalDefinitions"->True}
];

bifurMap=plotPoints[points,imageWidth,imageHeight]*N[((255*1.5)/seriesDepth)];
Dimensions[bifurMap]

image=ImageResize[Image[bifurMap],{2560,1080},Resampling->{"Gaussian",16,1/4}];
ImageResize[image,Large,Resampling->"Gaussian"]

Export[FileNameJoin[{NotebookDirectory[],"logisticmap.png"}], image];