くぴんのブログ

pythonの勉強としてやってみた。

↓条件
1. 初期条件→総人口100万人(男50万人、女50万人)、0~100歳までrandint(0,101)で一様分布。
2. 期間は100年間
3. 途中で死亡することはなく、すべて100歳まで生き、101歳になったら消える
4. 適齢期(18~40)の男か女、少ない方の数に割合(0.1など)をかけた数が毎年出生する
5. 男児の出生率を0.51

↓定常シナリオ (4の割合を0.1)


↓減少シナリオ (4の割合を0.02)




↓早婚 (15~35)、晩婚 (25~45)の影響を増加シナリオ (4の割合を0.2)で確認
↓早婚 (15~35)

↓晩婚 (25~45)






このモデルでは適齢期を約20年であり一様分布なため適齢期人口は20万人。カップル数は10万。ここから0.1の割合で出生すると1万人となって、ほぼ定常状態となる。

当たり前だけれど、増加シナリオでは、早婚であれば、再生産のサイクルが早く回るため、晩婚に比べて人口が早く増える。

このモデルでは、個体識別せずに、ざっくりと適齢期の半分をカップルとし、そこからある割合で再生産されるとしている。人口を増やすという観点からは、やはり結婚という形態にとらわれずに、適齢期の人たちが、毎年、前の年よりも多く生んでくれると効果がありそうだ。

可視化は以下のサイトを参考にした。
https://note.com/michell72/n/nf4183efdc5da

pngからgifを作成するには、以下のコマンドをsubprocessで実行した。環境はUbuntu 22.04.2 LTS。
cmd="convert -delay 30 -loop 0 graph0*png graph_.gif"
subprocess.call(cmd.split())

↓ソースコード
-----
"""
population simulation
"""
import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
import time
import subprocess

sta=time.time()

timespan=100
ratio=0.1
male_ratio=0.51
# children are born from couple between lower and upper years old.
lower=18
upper=40
lifespan=100
ini_pop=500000

print("timespan=",timespan,"ratio=",ratio,"male_ratio=",male_ratio,
        "lower=",lower,"upper=",upper,"lifespan=",lifespan,"ini_pop=",ini_pop)


def plot(male, female):
    fig, ax = plt.subplots(ncols=2, figsize=(6,8))

    male1=pd.Series([len(male[male==i]) for i in range(lifespan+1)])
    female1=pd.Series([len(female[female==i]) for i in range(lifespan+1)])

    # male
    ax[0].barh(male1.index, male1, color='b', height=0.7, label='male', alpha = 1)
    ax[0].yaxis.tick_right() # x reverse
    ax[0].set_yticklabels([]) #hide ticks
    ax[0].set_xlim([male1.max()*1.1, 0])# x reverse
    ax[0].set_title('MALE',fontsize=18)
    ax[0].xaxis.grid(True, which = 'major', linestyle = '--', color = '#CFCFCF')
    #plt.setp(ax[0].get_xticklabels(), fontsize=14)

    # female
    ax[1].barh(female1.index, female1, color='r', height=0.7, label='female', alpha = 1)
    ax[1].set_xlim([0,female1.max()*1.1])
    ax[1].set_title('FEMALE',fontsize=18)
    ax[1].xaxis.grid(True, which = 'major', linestyle = '--', color = '#CFCFCF')

    fig.suptitle("year="+str(year)+", population="+str(len(male)+len(female)), size=18)

    plt.savefig('graph'+str(year).rjust(4,"0")+'.png')
    #plt.show()
    plt.close(fig)

#initial male and female for 500,000
male=pd.Series(np.random.randint(0,lifespan+1,ini_pop))
female=pd.Series(np.random.randint(0,lifespan+1,ini_pop))

new_b=0

for year in range(timespan+1):
    plot(male, female)
    print("year=", year, "new_b=", new_b, "male= ",len(male), "female= ",
            len(female), "population= ", len(male)+len(female))
    # children are born from couple between lower and upper years old.
    young_m=len(male[lower<male][male<upper])
    young_f=len(female[lower<female][female<upper])
    new_b=int(ratio*min(young_m, young_f))
    n_male=int(male_ratio*new_b)
    # aging
    male=male+1
    female=female+1
    # children are added
    male=pd.concat([male, pd.Series(np.zeros(n_male))])
    female=pd.concat([female, pd.Series(np.zeros(new_b-n_male))])
    male=male.reset_index(drop=True)
    female=female.reset_index(drop=True)
    male=male[male<lifespan]
    female=female[female<lifespan]
plot(male, female)

cmd="convert -delay 30 -loop 0 graph0*png graph_.gif"
subprocess.call(cmd.split())

fin=time.time()
print("elapsed time [sec]=", fin-sta)
-----
何かの参考になれば幸いです。



にほんブログ村