## ----ins10--------------------------------------------------------------------
get_payout <- function(y, y_z, trigger, price){
   pmax( (trigger * y_z) - y,  0) * price
}


## ----ins11--------------------------------------------------------------------
yields <- seq(500, 2500, 100)
mean_yield <- mean(yields)
yields
mean_yield


## ----ins12--------------------------------------------------------------------
price <- 0.23


## ----ins13--------------------------------------------------------------------
payout <- get_payout(yields, mean_yield, trigger=0.9, price=price)
plot(yields, payout, type="l", xlab="yield", ylab="payout", las=1)


## ----ins112-------------------------------------------------------------------
rice_income <- yields * price
total_income <- rice_income + payout
plot(yields, total_income, type="l", xlab="yield", ylab="income", las=1)


## ----ins121, message=FALSE----------------------------------------------------
library(agrodata)
agro::ce_income(rice_income, rho=1.5)


## ----ins212, message=FALSE----------------------------------------------------
agro::ce_income(total_income, rho=1.5)


## ----ins20--------------------------------------------------------------------
rice_income <- yields * price
insurance_income <- payout - 50
total_income <- rice_income + insurance_income
agro::ce_income(rice_income, rho=1.5)
agro::ce_income(total_income, rho=1.5)


## ----ins30--------------------------------------------------------------------
plot(yields, rice_income, type="l", xlab="yield", ylab="income", las=1, ylim=c(-60,600), lwd=2)
lines(yields, insurance_income, col="blue", lty=2, lwd=2)
lines(yields, total_income, col="red", lwd=2)
legend("topleft", c("rice", "insurance", "total"), col=c("black", "blue", "red"), lty=c(1,2,1), lwd=2)


## ----ins40--------------------------------------------------------------------
rho <- seq(0,3,.1)
ce_rice <- sapply(rho, function(r) agro::ce_income(rice_income, rho=r))
ce_tot <- sapply(rho, function(r) agro::ce_income(total_income, rho=r))
plot(rho, ce_rice, col="red", type="l", las=1, ylab="CE", lwd=2)
lines(rho, ce_tot, col="blue", lwd=2)
legend("topright", c("no inusrance", "insurance"), col=c("red", "blue"), lty=1, lwd=2 )

# certainty equivalents are the same at
i <- which.min(abs(ce_tot - ce_rice))
rho[i]
abline(v=rho[i], lty=3, col="gray")


## ----ins2---------------------------------------------------------------------
z <- readRDS("hh_rice_yield.rds")
head(z)


## ----ins3---------------------------------------------------------------------
z$payout <- get_payout(z$y_zt, z$y_z, 0.9, price)
head(z)


## ----ins4---------------------------------------------------------------------
pay_ha_year <- tapply(z$payout, z$year, mean)
pay_ha_year
afp <- mean(pay_ha_year)
afp


## ----ins444-------------------------------------------------------------------
premium <- round(afp, 1) * 1.2
premium


## ----ins440-------------------------------------------------------------------
zz <- aggregate(z[, "payout", drop=FALSE], z[, "zone", drop=FALSE], mean)
zz$premium <- zz$payout * 1.2
zz


## ----ins50, out.width = '100%'------------------------------------------------
z$income <- z$y * price
z$income_with_ins <- z$income + z$payout - premium
head(z)


## ----ins06, out.width = '100%'------------------------------------------------
hh <- aggregate(z[, c("income", "income_with_ins")], z[,"fid", drop=FALSE], function(i) agro::ce_income(i, 1.5))

hh$benefit <- (hh$income_with_ins - hh$income)
hh$rel_benefit <- hh$benefit/ hh$income
head(hh)

mean(hh$benefit)
quantile(hh$benefit)


## ----ins066, out.width = '100%'-----------------------------------------------
zz <- aggregate(z[, c("income", "income_with_ins")], z[,"zone", drop=FALSE], function(i) agro::ce_income(i, 1.5))

zz$benefit <- (zz$income_with_ins - zz$income)
zz$rel_benefit <- zz$benefit/ zz$income
head(zz)

mean(zz$benefit)
quantile(zz$benefit)
barplot(sort(zz$rel_benefit), ylab="relative benefit", xlab="zone", las=1)