# R script used in https://remcat.hatenadiary.jp/entry/20210920/workerpop
# URL: http://tsigeto.info/maikin/maikin-monthly2.r.txt
# 2021/09/14 - 2021/09/19
# Created by TANAKA Sigeto <http://tsigeto.info/maikin/>
# Run http://tsigeto.info/maikin/maikin-monthly.r.txt before using it.

for( s in names(data) ) {
	n <- nrow( data[[s]] )
	data[[s]]$worker.prev <- c( NA, data[[s]] [ -n, "worker" ] )
	data[[s]]$worker.next <- c(     data[[s]] [ -1, "worker" ] , NA )
	data[[s]]$worker.inc <- log( data[[s]]$worker / data[[s]]$worker.prev )
	data[[s]] [ data[[s]]$yyyymm %in% reset.date1 , "worker.inc" ] <- NA
}

worker.cum <- function( r , start ) {
	r[ is.na(r$worker.inc), "worker.inc" ] <- 0
	r$cum <- cumsum( r$worker.inc )
	r$worker.cum <- exp(r$cum) * start
	r
}
cum0 <- lapply(
	data,
	function(d) { 
		start <- d[1, "worker"]
		worker.cum( subset( d, yyyymm %% 1 == 0  ) , start )
	}
)
cum5 <- lapply(
	data,
	function(d) { 
		start <- d[1, "worker"]
		worker.cum( subset( d, yyyymm %% 1 == 0.5  ) , start )
	}
)

temp <- sapply (
	cum0[[1]]$yyyymm ,
	function(i) {
		j <- as.character(i)
		k <- as.character(i-0.5)
		w <- rownames(worker.pop)
		if (     j %in% w ) { r <- worker.pop[ j, ] }
		else if( k %in% w ) { r <- worker.pop[ k, ] }
		else { r <- rep( NA, ncol(worker.pop) ) }
		r
	}	
)
worker.pop.yyyymm <- t(temp)
colnames(worker.pop.yyyymm) <- colnames(worker.pop)
rownames(worker.pop.yyyymm) <- cum0[[1]]$yyyymm

worker.predicted <- list()
for( s in names(cum5) ) {
	size <- cum5[[s]] [1, "size"]
	worker.predicted[[s]] <- subset( cum0[[s]] , TRUE,  select="yyyymm" )
	worker.predicted[[s]]$e1.e0 <- cum0[[s]] [ , "worker.cum"] 
	worker.predicted[[s]]$e0.e1 <- cum5[[s]] [ , "worker.cum"] 
	worker.predicted[[s]]$census <- worker.pop.yyyymm[ , as.character(size) ]
}

renew.bm <- function( d , e ) {
	prev <- c( NA, d[ -nrow(d) , e ] )
	d$gap <- d[ , e ] / prev

	d$worker.bm <- NA
	d[1, "worker.bm"] <- d[ 1, e ]
	d[1, "gap"] <- 1

	for( i in 2:nrow(d) ) {
		gap <-  d[ i, "gap" ]
		mon <- d[i,"yyyymm"]
		if( 200901 == mon ) { gap <- g2009 }
		else if(201201==mon){ gap <- g2012 }
		else if(201801==mon){ gap <- g2018 }
		d[ i, "worker.bm" ] <- d[ i-1, "worker.bm"] * gap
		d[ i, "gap" ] <- gap

		g <- d[ i, "census" ] / d[ i, "worker.bm" ]
		if( 200610 == mon ) { g2009 <- g }
		else if(200907==mon){ g2012 <- g }
		else if(201407==mon){ g2018 <- g } 
	}
	d
}

with.bm.e1.e0 <- lapply( worker.predicted , renew.bm , "e1.e0" )
with.bm.e0.e1 <- lapply( worker.predicted , renew.bm , "e0.e1" )

with.bm <- with.bm.e1.e0
for( s in names(with.bm)  ) {
	with.bm[[s]] <- subset( with.bm[[s]] , TRUE, select="yyyymm" )
	with.bm[[s]]$e1.e0 <- with.bm.e1.e0[[s]]$worker.bm
	with.bm[[s]]$e0.e1 <- with.bm.e0.e1[[s]]$worker.bm
	with.bm[[s]]$census <-with.bm.e0.e1[[s]]$census
}

inc.prev12mon <- function(v) {
	p <- c( rep(NA,12) , v[ 1:( length(v)-12 ) ] )
	v / p
}

inc12mon <- lapply( 
	worker.predicted,
	function(d){
		r1 <- inc.prev12mon( d$e1.e0 )
		r2 <- inc.prev12mon( d$e0.e1 )
		r <- cbind( r1 , r2 )
		colnames( r ) <- c( "e1.e0", "e0,e1" ) 
		rownames( r ) <- d$yyyymm
		data.frame(r)
	}
)

inc12mon.real.e1 <- lapply (
	data ,
	function(d) {
		start <- d[1, "worker"]
		wc <- worker.cum( d , start )
		wc <- subset( wc, yyyymm %% 1 == 0.5 , select="worker.cum" )
		inc.prev12mon(wc$worker.cum)
	}
)

for ( i in 1:length(inc12mon) ) {
	inc12mon[[i]] $ e1.real <- inc12mon.real.e1[[i]]
}

interval.mon <- c(
	33,	# 200401-200610
	33,	# 200610-200907
	60,	# 200907-201407
	23,	# 201407-201606
	36,	# 201606-201906
	23	# 201906-202105
)

worker.pop.inc <- apply(
	worker.pop,
	2,
	function(v) {
		p <- c( NA, v )
		n <- c( v, NA )
		n/p
	}
)

worker.pop.inc.yearly <- exp( log(worker.pop.inc) * 12 / interval.mon )

yyyymm <- cum0[[1]]$yyyymm
temp <- matrix( nrow=length(yyyymm), ncol=ncol(worker.pop.inc) )
rownames(temp)<-yyyymm
colnames(temp)<-names(cum0)
for( end in rev( rownames(worker.pop.inc) )[-1] ) {
	for( size in 1:ncol(worker.pop.inc)) {
		temp[ yyyymm<=end, size ] <- worker.pop.inc.yearly[end,size]
	}
}
log12mon.census <- list()
for( s in names(inc12mon) ){
	log12mon.census[[s]] <- log( inc12mon[[s]] )
	log12mon.census[[s]] $ census <- log( temp[ , s] )
}