Here we will investigate the nontrivial equivlibrium solutions to the α model for specific parameter values.

(* Define the birth functions *)
RickerB[x_,y_]:=Exp[-(x+y)];

(* Define the parameters *)
(* Both Ponds *)
pS=0.4;
pI=0.3;
mS=0.3;
mI=0.2;
r=0.15;
(* Pond 1 *)
w1=1;
(* Pond 2 *)
w2=1;

(* Define the itteration functions *)

Sus1[S1_,I1_,S2_,I2_,b_,α_]:=b*RickerB[S1, I1]*S1 + pS*Exp[-w1*I1]*S1+pI*r*I1+α*mS*Exp[-w2*I2]*S2+α*mI*r*I2;
Inf1[S1_,I1_,S2_,I2_,b_,α_]:=b*RickerB[S1, I1]*I1 + pS*(1-Exp[-w1*I1])*S1 + pI*(1-r)*I1+α*mS*(1-Exp[-w2*I2])*S2+α*mI*(1-r)*I2;
Sus2[S1_,I1_,S2_,I2_,b_,α_]:=α*b*RickerB[S2, I2]*S2 + α*pS*Exp[-w2*I2]*S2+α*pI*r*I2+mS*Exp[-w1*I1]*S1+mI*r*I1;
Inf2[S1_,I1_,S2_,I2_,b_,α_]:=α*b*RickerB[S2, I2]*I2 + α*pS*(1-Exp[-w2*I2])*S2 + α*pI*(1-r)*I2+mS*(1-Exp[-w1*I1])*S1+mI*(1-r)*I1;

EqSols[x0_,y0_,c_,n_]:=
Module[{b=x0,α=y0,i,S1i,I1i,S2i,I2i,S1n,I1n,S2n,I2n,S1a,I1a,S2a,I2a},
S1i=b*(1.5-Random[]);
I1i=b*(1.5-Random[]);
S2i=b*(1.5-Random[]);
I2i=b*(1.5-Random[]);
For[ i=1,i<c,i++,
S1n=Sus1[S1i,I1i,S2i,I2i,b,α];
I1n=Inf1[S1i,I1i,S2i,I2i,b,α];
S2n=Sus2[S1i,I1i,S2i,I2i,b,α];
I2n=Inf2[S1i,I1i,S2i,I2i,b,α];
S1i=S1n;
I1i=I1n;
S2i=S2n;
I2i=I2n;
];
S1a=0.0;
I1a=0.0;
S2a=0.0;
I2a=0.0;
For[ i=1,i<n,i++,
S1n=Sus1[S1i,I1i,S2i,I2i,b,α];
I1n=Inf1[S1i,I1i,S2i,I2i,b,α];
S2n=Sus2[S1i,I1i,S2i,I2i,b,α];
I2n=Inf2[S1i,I1i,S2i,I2i,b,α];
S1i=S1n;
I1i=I1n;
S2i=S2n;
I2i=I2n;
S1a=S1a+S1i/n;
I1a=I1a+I1i/n;
S2a=S2a+S2i/n;
I2a=I2a+I2i/n;
];
Outv=FindRoot[{b*RickerB[S1, I1]*S1 + pS*Exp[-w1*I1]*S1+pI*r*I1+α*mS*Exp[-w2*I2]*S2+α*mI*r*I2==S1,b*RickerB[S1, I1]*I1 + pS*(1-Exp[-w1*I1])*S1 + pI*(1-r)*I1+α*mS*(1-Exp[-w2*I2])*S2+α*mI*(1-r)*I2==I1,
α*b*RickerB[S2, I2]*S2 + α*pS*Exp[-w2*I2]*S2+α*pI*r*I2+mS*Exp[-w1*I1]*S1+mI*r*I1==S2,α*b*RickerB[S2, I2]*I2 + α*pS*(1-Exp[-w2*I2])*S2 + α*pI*(1-r)*I2+mS*(1-Exp[-w1*I1])*S1+mI*(1-r)*I1==I2},
{{S1,S1a},{I1,I1a},{S2,S2a},{I2,I2a}}] ;
Outv 
]

EqSols[5,1.2,100,2]

beginα = 0;
endα = Min[(1-mS)/pS,(1-mI)/pI,(1-pS)/mS,(1-pI)/mI];

JacMat[f_,x_]:=Outer[D,f,x]

JM[b_,α_,c_,n_]:=JacMat[{Sus1[S1,I1,S2,I2,b,α],Inf1[S1,I1,S2,I2,b,α],Sus2[S1,I1,S2,I2,b,α],Inf2[S1,I1,S2,I2,b,α]},{S1,I1,S2,I2}]/.EqSols[b,α,c,n]

MatrixForm[JM[5,1.2,100,2]]

Eigenvalues[JM[60,0.8,100,2]]

S1EqPlot=Plot[S1/.EqSols[5,α,100,2],{α,beginα,endα} ,AxesLabel->{"α", "(Subscript[S, 1])"},AxesOrigin->{0,0}];
I1EqPlot=Plot[I1/.EqSols[5,α,100,2],{α,beginα,endα},AxesLabel->{"α", "(Subscript[I, 1])"},AxesOrigin->{0,0}];
S2EqPlot=Plot[S2/.EqSols[5,α,100,2],{α,beginα,endα},AxesLabel->{"α", "(Subscript[S, 2])"},AxesOrigin->{0,0}];
I2EqPlot=Plot[I2/.EqSols[5,α,100,2],{α,beginα,endα},AxesLabel->{"α", "(Subscript[I, 2])"},AxesOrigin->{0,0}];
GraphicsGrid[{{S1EqPlot, I1EqPlot},{S2EqPlot, I2EqPlot}},Spacings->{Scaled[0.25],Scaled[0.25]}]

Export["T:/white/research/ICTCM/NumEqB5S1_1.pdf",S1EqPlot]
Export["T:/white/research/ICTCM/NumEqB5I1_1.pdf",I1EqPlot]
Export["T:/white/research/ICTCM/NumEqB5S2_1.pdf",S2EqPlot]
Export["T:/white/research/ICTCM/NumEqB5I2_1.pdf",I2EqPlot]

LamPlot=Plot[Abs[Eigenvalues[JM[5,a,100,2]]],{a,beginα,endα} ,AxesLabel->{"α","|λ|"},AxesOrigin->{0,0}]

Export["T:/white/research/ICTCM/NumLamB5.pdf",LamPlot]

S1EqPlot=Plot[S1/.EqSols[20,α,200,20],{α,beginα,endα} ,AxesLabel->{"α", "(Subscript[S, 1])"},AxesOrigin->{0,0}];
I1EqPlot=Plot[I1/.EqSols[20,α,200,20],{α,beginα,endα},AxesLabel->{"α", "(Subscript[I, 1])"},AxesOrigin->{0,0}];
S2EqPlot=Plot[S2/.EqSols[20,α,200,20],{α,beginα,endα},AxesLabel->{"α", "(Subscript[S, 2])"},AxesOrigin->{0,0}];
I2EqPlot=Plot[I2/.EqSols[20,α,200,20],{α,beginα,endα},AxesLabel->{"α", "(Subscript[I, 2])"},AxesOrigin->{0,0}];
GraphicsGrid[{{S1EqPlot, I1EqPlot},{S2EqPlot, I2EqPlot}},Spacings->{Scaled[0.25],Scaled[0.25]}]

Export["T:/white/research/ICTCM/NumEqB20S1_1.pdf",S1EqPlot]
Export["T:/white/research/ICTCM/NumEqB20I1_1.pdf",I1EqPlot]
Export["T:/white/research/ICTCM/NumEqB20S2_1.pdf",S2EqPlot]
Export["T:/white/research/ICTCM/NumEqB20I2_1.pdf",I2EqPlot]

LamPlot=Plot[Abs[Eigenvalues[JM[20,a,200,20]]],{a,beginα,endα} ,AxesLabel->{"α","|λ|"},AxesOrigin->{0,0}]

Export["T:/white/research/ICTCM/NumLamB20.pdf",LamPlot]

S1EqPlot=Plot[S1/.EqSols[60,α,200,20],{α,beginα,endα} ,AxesLabel->{"α", "(Subscript[S, 1])"},AxesOrigin->{0,0}];
I1EqPlot=Plot[I1/.EqSols[60,α,200,20],{α,beginα,endα},AxesLabel->{"α", "(Subscript[I, 1])"},AxesOrigin->{0,0}];
S2EqPlot=Plot[S2/.EqSols[60,α,200,20],{α,beginα,endα},AxesLabel->{"α", "(Subscript[S, 2])"},AxesOrigin->{0,0}];
I2EqPlot=Plot[I2/.EqSols[60,α,200,20],{α,beginα,endα},AxesLabel->{"α", "(Subscript[I, 2])"},AxesOrigin->{0,0}];
GraphicsGrid[{{S1EqPlot, I1EqPlot},{S2EqPlot, I2EqPlot}},Spacings->{Scaled[0.25],Scaled[0.25]}]

Export["T:/white/research/ICTCM/NumEqB60S1_1.pdf",S1EqPlot]
Export["T:/white/research/ICTCM/NumEqB60I1_1.pdf",I1EqPlot]
Export["T:/white/research/ICTCM/NumEqB60S2_1.pdf",S2EqPlot]
Export["T:/white/research/ICTCM/NumEqB60I2_1.pdf",I2EqPlot]

LamPlot=Plot[Abs[Eigenvalues[JM[60,a,200,20]]],{a,beginα,endα} ,AxesLabel->{"α","|λ|"},AxesOrigin->{0,0}]

Export["T:/white/research/ICTCM/NumLamB60.pdf",LamPlot]