齢構造を持つモデル

Leslie matrix models

ハイイロリスのデータ

The population consists of 7 age classes.

leslie = {{f1,f2,f3,f4,f5,f6,f7},
		  {P1,0,0,0,0,0,0},
		  {0,P2,0,0,0,0,0},
		  {0,0,P3,0,0,0,0},
		  {0,0,0,P4,0,0,0},
		  {0,0,0,0,P5,0,0},
		  {0,0,0,0,0,P6,0}}

MatrixForm[leslie]

Data of grey squirrels in North Carolina

para = {f1->0.32, f2->0.57, f3->0.57, f4->0.57,
	    f5->0.57, f6->0.57,f7->0.57,
	    P1->0.46, P2->0.77, P3->0.65, P4->0.67,
	    P5->0.64, P6->0.88};

leslie2 = leslie/.para

初期分布1: 最若齢個体から開始

init = {0.1,0,0,0,0,0,0};

leslie2.init

leslie2.leslie2.init

MatrixPower[leslie2,2].init

MatrixPower[leslie2,3].init

MatrixPower[leslie2,4].init

dist[t_]:=MatrixPower[leslie2,t].init

dist[10]

dist[20]

dist[100]

data = Table[dist[t],{t,0,20}]

BarChart[data[[3]],PlotRange->{-0.01,0.1},ChartLabels->{"1","2","3","4","5","6","7"}]

?BarChart

glist={};
Do[
g=BarChart[data[[i]],PlotRange->{-0.01,0.1},ChartLabels->{"1","2","3","4","5","6","7"}];
AppendTo[glist,g],{i,1,Length[data]}
]

glist

ListAnimate[glist]

Export["AgeDistDynamics.mov",glist]

data = Table[dist[t],{t,0,100}];

totalDensity = Apply[Plus,Transpose[data]]

ListPlot[totalDensity,Joined->True,PlotRange->All,PlotStyle->Thick,AxesLabel->{"t","PopSize Total"}]

ListLogPlot[totalDensity,Joined->True,PlotStyle->Thick,AxesLabel->{"t","PopSize Total"}]

Fit[ Log[totalDensity],{1,t},t]

Exp[0.0377887]

初期分布2: 最老齢個体から開始

init = {0,0,0,0,0,0,0.1};

leslie2.init

leslie2.leslie2.init

MatrixPower[leslie2,2].init

MatrixPower[leslie2,3].init

MatrixPower[leslie2,4].init

dist[t_]:=MatrixPower[leslie2,t].init

dist[10]

dist[20]

dist[100]

data2 = Table[dist[t],{t,0,20}]

BarChart[data2[[1]],PlotRange->{-0.01,0.1},ChartLabels->{"1","2","3","4","5","6","7"}]

glist2={};
Do[
g=BarChart[data2[[i]],PlotRange->{-0.01,0.1},ChartLabels->{"1","2","3","4","5","6","7"}];
AppendTo[glist2,g],{i,1,Length[data2]}
]

glist2

ListAnimate[glist2]

Export["AgeDistDynamics2.mov",glist2]

data2 = Table[dist[t],{t,0,100}];

totalDensity2 = Apply[Plus,Transpose[data2]]

ListPlot[totalDensity2,Joined->True,PlotRange->All,PlotStyle->Thick,AxesLabel->{"t","PopSize Total"}]

ListLogPlot[totalDensity2,Joined->True,PlotStyle->Thick,AxesLabel->{"t","PopSize Total"}]

ListLogPlot[{totalDensity,totalDensity2},Joined->True,PlotStyle->Thick,AxesLabel->{"t","PopSize Total"}]

Fit[ Log[totalDensity],{1,t},t]

Exp[0.0377887]

行列 leslie2 の固有値を求める

Eigenvalues[leslie2]//Chop

行列 leslie2 の固有ベクトルを求める

Eigenvectors[leslie2]//Chop

ratios = %[[1]]

ratios/Apply[Plus,ratios]

Last[data]

Last[data]/Apply[Plus, Last[data] ]

fertility = f1 + f2 P1 + f3 P1 P2 + f4 P1 P2 P3 + f5 P1 P2 P3 P4 +
	f6 P1 P2 P3 P4 P5 + f7 P1 P2 P3 P4 P5 P6

fertility/.para

周期行列

mat={{0,0,6},{1/2,0,0},{0,1/3,0}}

eigens=Eigenvalues[mat]

gEigens=ComplexListPlot[eigens,PlotStyle->Automatic,PlotRange->{{-1,1},{-1,1}},AxesLabel->{"Re","Im"}]

gCircle=ParametricPlot[{Cos[theta],Sin[theta]},{theta,0,2Pi}]

Show[gEigens,gCircle]

Eigenvectors[mat]

init={1,0,0}

dist[t_]:=MatrixPower[mat,t].init

dist[1]

data = Table[dist[t],{t,0,20}]

glist={};
Do[
g=BarChart[data[[i]],PlotRange->{-0.01,1.2},ChartLabels->{"1","2","3","4","5","6","7"}];
AppendTo[glist,g],{i,1,Length[data]}
]

glist[[1]]

glist[[2]]

glist[[3]]

ListAnimate[glist]

totalDensity = Apply[Plus,Transpose[data]]

ListPlot[totalDensity,Joined->True,PlotRange->All,PlotStyle->Thick,AxesLabel->{"t","PopSize Total"}]

密度依存する行列

mat={{sigma1(1-gamma),p0 Exp[-b(n1+n2)]},{sigma1 gamma, sigma2}}

para={b->1,sigma1->0.5,sigma2->0.1,gamma->0.1}

matN=mat/.para/.p0->300

nextStatus[{n1x_,n2x_}]:=matN.{n1,n2}/.{n1->n1x,n2->n2x}

nextStatus[{0.1,0.2}]

seq=NestList[nextStatus,{0.1,0.2},50]

totalDensity=Apply[Plus,Transpose[seq]]

ListPlot[totalDensity,Joined->True,PlotStyle->Thick,AxesLabel->{"t","n1+n2"},PlotRange->All]

Mathematicaの行列計算

2 x2行列の定義

mat={{1,2},{0,-3}}

MatrixForm[mat]

行列同士のかけ算

mat.mat

MatrixPower[mat,2]

MatrixPower[mat,5]

固有値

Eigenvalues[mat]

固有値と固有ベクトル

Eigensystem[mat]

?Eigensystem

mat2={{1,2,3},{0,9,4},{-1,0,4}}

MatrixForm[mat2]

Eigenvalues[mat2]//N