擬似健康保険データを作ってみる。 会員全体のデータを作るために、出生と死亡のシミュレーションを行う。

まずは出生。 年齢と性別に対して何人生まれたか…というデータを、擬似的に生成する。

IPSSが公開している、年次出生数のデータを利用する。

library(readr)
library(dplyr)
url = 'http://www.ipss.go.jp/p-toukei/JMD/00/STATS/Births.txt'
dat = read.table(url, skip=2, header = TRUE)
dat %>% head
##   Year  Female    Male   Total
## 1 1947 1301806 1376986 2678792
## 2 1948 1303060 1378564 2681624
## 3 1949 1316630 1380008 2696638
## 4 1950 1134396 1203111 2337507
## 5 1951 1043048 1094641 2137689
## 6 1952  977101 1028061 2005162

ある年に、FemaleとMaleが何人生まれたかのデータになっている。

性別列を作って縦に持っておく。

library(tidyr)
library(magrittr)
dat %>%  
 pivot_longer(cols = c("Male", "Female"), names_to = "Sex", values_to = "Life") %>% 
 mutate(Sex = if_else(Sex == "Female", "F", "M")) -> dat
dat
## # A tibble: 140 x 4
##     Year   Total Sex      Life
##    <int>   <int> <chr>   <int>
##  1  1947 2678792 M     1376986
##  2  1947 2678792 F     1301806
##  3  1948 2681624 M     1378564
##  4  1948 2681624 F     1303060
##  5  1949 2696638 M     1380008
##  6  1949 2696638 F     1316630
##  7  1950 2337507 M     1203111
##  8  1950 2337507 F     1134396
##  9  1951 2137689 M     1094641
## 10  1951 2137689 F     1043048
## # … with 130 more rows

乱数でシミュレーションをするために、累積比率を出しておく

dat %<>%
 arrange(Sex, Year)
dat %<>%
 mutate(ratio = Life / sum(Life)) %>% 
 mutate(cum_sum = cumsum(ratio))

dat %>% head
## # A tibble: 6 x 6
##    Year   Total Sex      Life   ratio cum_sum
##   <int>   <int> <chr>   <int>   <dbl>   <dbl>
## 1  1947 2678792 F     1301806 0.0121   0.0121
## 2  1948 2681624 F     1303060 0.0122   0.0243
## 3  1949 2696638 F     1316630 0.0123   0.0366
## 4  1950 2337507 F     1134396 0.0106   0.0472
## 5  1951 2137689 F     1043048 0.00973  0.0569
## 6  1952 2005162 F      977101 0.00912  0.0660
dat %>% 
 select(-Total) %>% 
 write.csv("ipss_birth.csv", row.names=FALSE, quote = FALSE)

累積確率分布を作る。これにより、一様乱数(y軸:0-1)から、対応する行番号が得られて、 行番号から、(出生年,性別)を得ることができる。

欲しい人数分だけ一様乱数を取れば、出生のデータが得られる。

library(ggplot2)
dat = read_csv("ipss_birth.csv")
## Parsed with column specification:
## cols(
##   Year = col_double(),
##   Sex = col_character(),
##   Life = col_double(),
##   ratio = col_double(),
##   cum_sum = col_double()
## )
dat %>% 
 mutate(x_axis = paste(Year, Sex, sep="")) %>% 
 mutate(row_num = 1:nrow(.)) %>% 
 ggplot() + 
 geom_line(aes(x = row_num, y = cum_sum), stat = "identity") + 
 annotate("segment",x=1,xend=90,y=0.7,yend=0.7,colour="blue",
             size=1, arrow=arrow()) + 
 annotate("segment",x=92,xend=92,y=0.65,yend=0.0,colour="blue",
             size=1,arrow=arrow()) + 
 annotate("text", x=20, y=0.75, parse=TRUE, label="'Random Number 0.7'") +
 annotate("text", x=120, y=0.2, parse=TRUE, label="'Rownum: 95 -> 1969-Female'") +
 ggtitle("Cumulative ratio of Population Birth")

これを利用して、出生のシミュレーションができる。

次は、死亡率を算出して、死亡のシミュレーションを行う。