pacman::p_load(
tidyverse,
kableExtra,
latex2exp
)4 Comparison: Risk sharing with limited commitment but without storage
Finally, I compare the model with a risk-sharing model under limited commitment but without storage, as in Kocherlakota (1996) and Ligon, Thomas, and Worrall (2002). I use the simulation code I wrote on this webstie. After I show that, in contrast to the model with saving, the model here gives consumption whcih takes a finite number of values in limit, I demonstrate that the welfare actually can decrease due to the existence of storage.
4.1 Code
4.1.1 Utility functions
calculateUtility <- function(cons, sigma) {
if (sigma != 1) {
utility = (cons^(1 - sigma) - 1) / (1 - sigma)
} else if (sigma == 1) {
utility = log(cons)
}
return(utility)
}
calculateMarginalUtility <- function(cons, sigma) cons^(- sigma)4.1.2 Consumption function
calculateHH1Consumption <- function(
aggregateResources,
relativeParetoWeight,
sigma
) {
aggregateResources / (1 + (relativeParetoWeight^(1 / sigma)))
}4.1.3 Value under autarky
calculateAutarkyValue <- function(
incomeGridPoints,
sigma,
delta,
punishment,
incomeTransitionMatrix
) {
autarkyValue <- numeric(length = length(incomeGridPoints))
i <- 1
diff <- 1
while (diff > 1e-12) {
autarkyValueNew <- (
calculateUtility(incomeGridPoints * (1 - punishment), sigma)
+ delta * incomeTransitionMatrix %*% autarkyValue
)
diff <- max(abs(autarkyValueNew - autarkyValue))
autarkyValue <- autarkyValueNew
i <- i + 1
}
return(autarkyValue)
}getRelativeParetoWeightsGridPoints <- function(
sigma,
punishment,
householdIncomeGridPoints,
villageIncomeGridPoints,
numRelativeParetoWeights
) {
minRelativeParetoWeights <- (
calculateMarginalUtility(max(villageIncomeGridPoints), sigma)
/ calculateMarginalUtility(min(householdIncomeGridPoints * (1 - punishment)), sigma)
)
maxRelativeParetoWeights <- (
calculateMarginalUtility(min(villageIncomeGridPoints * (1 - punishment)), sigma)
/ calculateMarginalUtility(max(householdIncomeGridPoints), sigma)
)
relativeParetoWeightsGridPoints <- exp(
seq(
log(minRelativeParetoWeights),
log(maxRelativeParetoWeights),
length.out = numRelativeParetoWeights)
)
return(relativeParetoWeightsGridPoints)
}calculateHouseholdConsumption <- function(
aggregateIncome,
relativeParetoWeight,
numHouseholds,
sigma
) {
aggregateIncome / (1 + (numHouseholds - 1) * (relativeParetoWeight^(- 1 / sigma)))
}calculateValueFullRiskSharing <- function(
incomeTransitionMatrix,
aggregateIncomeGridPoints,
delta,
sigma,
autarkyValueMatrix,
consumptionOnRelativeParetoWeightGrid,
numRelativeParetoWeights,
numHouseholds
) {
# Initial guess is expected utilities under autarky
householdValueFullRiskSharing <- outer(
autarkyValueMatrix[, 1], rep(1, numRelativeParetoWeights)
)
villageValueFullRiskSharing <- outer(
autarkyValueMatrix[, 2], rep(1, numRelativeParetoWeights)
)
iteration <- 1
diff <- 1
while (diff > 1e-10 & iteration < 500) {
householdValueFullRiskSharingNew <- (
calculateUtility(consumptionOnRelativeParetoWeightGrid, sigma)
+ delta * incomeTransitionMatrix %*% householdValueFullRiskSharing
)
villageValueFullRiskSharingNew <- (
calculateUtility(
(aggregateIncomeGridPoints - consumptionOnRelativeParetoWeightGrid) / (numHouseholds - 1),
sigma
)
+ delta * incomeTransitionMatrix %*% villageValueFullRiskSharing
)
diff <- max(
max(abs(householdValueFullRiskSharing - householdValueFullRiskSharingNew)),
max(abs(villageValueFullRiskSharing - villageValueFullRiskSharingNew))
)
householdValueFullRiskSharing <- householdValueFullRiskSharingNew
villageValueFullRiskSharing <- villageValueFullRiskSharingNew
iteration <- iteration + 1
}
return(list(
householdValueFullRiskSharing = householdValueFullRiskSharing,
villageValueFullRiskSharing = villageValueFullRiskSharing
))
}interpolateValueFunction <- function(
relativeParetoWeight,
relativeParetoWeightsGridPoints,
valueFunctionMatrix
) {
apply(
valueFunctionMatrix,
1,
function(x) {
approx(
relativeParetoWeightsGridPoints,
x,
relativeParetoWeight,
rule = 2
)$y
}
)
}
calculateDiffLCRiskSharingAndAutarky <- function(
relativeParetoWeight,
relativeParetoWeightsGridPoints,
delta,
sigma,
aggregateIncome,
householdValueLCRiskSharing,
villageValueLCRiskSharing,
incomeTransitionProbVec,
householdAutarkyValue,
villageAutarkyValue,
numHouseholds
) {
householdConsumption <- calculateHouseholdConsumption(
aggregateIncome,
relativeParetoWeight,
numHouseholds,
sigma
)
householdValueLCRiskSharingAtRelativeParetoWeight <- interpolateValueFunction(
relativeParetoWeight,
relativeParetoWeightsGridPoints,
householdValueLCRiskSharing
)
villageValueLCRiskSharingAtRelativeParetoWeight <- interpolateValueFunction(
relativeParetoWeight,
relativeParetoWeightsGridPoints,
villageValueLCRiskSharing
)
householdDiffLCRiskSharingAndAutarky <- (
calculateUtility(householdConsumption, sigma)
+ delta * incomeTransitionProbVec %*% householdValueLCRiskSharingAtRelativeParetoWeight
- householdAutarkyValue
) %>% as.numeric
villageDiffLCRiskSharingAndAutarky <- (
calculateUtility((aggregateIncome - householdConsumption) / (numHouseholds - 1), sigma)
+ delta * incomeTransitionProbVec %*% villageValueLCRiskSharingAtRelativeParetoWeight
- villageAutarkyValue
) %>% as.numeric
return(list(
householdDiffLCRiskSharingAndAutarky = householdDiffLCRiskSharingAndAutarky,
villageDiffLCRiskSharingAndAutarky = villageDiffLCRiskSharingAndAutarky
))
}
calculateValueLCRiskSharing <- function(
valueFullRiskSharing,
consumptionOnRelativeParetoWeightGrid,
aggregateIncomeGridPoints,
incomeTransitionMatrix,
autarkyValueMatrix,
relativeParetoWeightsGridPoints,
numRelativeParetoWeights,
delta,
sigma,
numIncomeStates,
numHouseholds,
iterationLimit,
diffLimit
) {
# Initial guess is expected utilities under full risk sharing
householdValueLCRiskSharing <- valueFullRiskSharing$householdValueFullRiskSharing
villageValueLCRiskSharing <- valueFullRiskSharing$villageValueFullRiskSharing
diff <- 1
iteration <- 1
while ((diff > diffLimit) && (iteration <= iterationLimit)) {
# First, ignore enforceability and just update the value functions
# using the values at the previous iteration
householdValueLCRiskSharingNew <- (
calculateUtility(consumptionOnRelativeParetoWeightGrid, sigma)
+ delta * incomeTransitionMatrix %*% householdValueLCRiskSharing
)
villageValueLCRiskSharingNew <- (
calculateUtility(
(aggregateIncomeGridPoints - consumptionOnRelativeParetoWeightGrid) / (numHouseholds - 1),
sigma
)
+ delta * incomeTransitionMatrix %*% villageValueLCRiskSharing
)
# Now check enforceability at each state
for (incomeStateIndex in seq(1, numIncomeStates)) {
householdAutarkyValue <- autarkyValueMatrix[incomeStateIndex, 1]
villageAutarkyValue <- autarkyValueMatrix[incomeStateIndex, 2]
if (any(householdValueLCRiskSharingNew[incomeStateIndex, ] <= householdAutarkyValue)) {
villageValueLCRiskSharingNew[
incomeStateIndex,
householdValueLCRiskSharingNew[incomeStateIndex, ] <= householdAutarkyValue
] <- villageValueLCRiskSharingNew[
incomeStateIndex,
householdValueLCRiskSharingNew[incomeStateIndex, ] <= householdAutarkyValue
] %>% min
householdValueLCRiskSharingNew[
incomeStateIndex,
householdValueLCRiskSharingNew[incomeStateIndex, ] <= householdAutarkyValue
] <- householdAutarkyValue
}
if (any(villageValueLCRiskSharingNew[incomeStateIndex, ] <= villageAutarkyValue)) {
householdValueLCRiskSharingNew[
incomeStateIndex,
villageValueLCRiskSharingNew[incomeStateIndex, ] <= villageAutarkyValue
] <- householdValueLCRiskSharingNew[
incomeStateIndex,
villageValueLCRiskSharingNew[incomeStateIndex, ] <= villageAutarkyValue
] %>% min
villageValueLCRiskSharingNew[
incomeStateIndex,
villageValueLCRiskSharingNew[incomeStateIndex, ] <= villageAutarkyValue
] <- villageAutarkyValue
}
}
diff <- max(
max(abs(householdValueLCRiskSharingNew - householdValueLCRiskSharing)),
max(abs(villageValueLCRiskSharingNew - villageValueLCRiskSharing))
)
householdValueLCRiskSharing <- householdValueLCRiskSharingNew
villageValueLCRiskSharing <- villageValueLCRiskSharingNew
iteration <- iteration + 1
}
relativeParetoWeightBounds <- matrix(NA, nrow = numIncomeStates, ncol = 2)
for (incomeStateIndex in seq(1, numIncomeStates)) {
aggregateIncome <- aggregateIncomeGridPoints[incomeStateIndex]
incomeTransitionProbVec <- incomeTransitionMatrix[incomeStateIndex,]
householdAutarkyValue <- autarkyValueMatrix[incomeStateIndex, 1]
villageAutarkyValue <- autarkyValueMatrix[incomeStateIndex, 2]
if (
calculateDiffLCRiskSharingAndAutarky(
min(relativeParetoWeightsGridPoints),
relativeParetoWeightsGridPoints,
delta,
sigma,
aggregateIncome,
householdValueLCRiskSharing,
villageValueLCRiskSharing,
incomeTransitionProbVec,
householdAutarkyValue,
villageAutarkyValue,
numHouseholds
)$householdDiffLCRiskSharingAndAutarky < 0) {
relativeParetoWeightLowerBound <- uniroot(
function(x) {calculateDiffLCRiskSharingAndAutarky(
x,
relativeParetoWeightsGridPoints,
delta,
sigma,
aggregateIncome,
householdValueLCRiskSharing,
villageValueLCRiskSharing,
incomeTransitionProbVec,
householdAutarkyValue,
villageAutarkyValue,
numHouseholds
)$householdDiffLCRiskSharingAndAutarky},
c(min(relativeParetoWeightsGridPoints), max(relativeParetoWeightsGridPoints)),
tol = 1e-10,
maxiter = 300
)$root
} else {
relativeParetoWeightLowerBound <- min(relativeParetoWeightsGridPoints)
}
if (
calculateDiffLCRiskSharingAndAutarky(
max(relativeParetoWeightsGridPoints),
relativeParetoWeightsGridPoints,
delta,
sigma,
aggregateIncome,
householdValueLCRiskSharing,
villageValueLCRiskSharing,
incomeTransitionProbVec,
householdAutarkyValue,
villageAutarkyValue,
numHouseholds
)$villageDiffLCRiskSharingAndAutarky < 0) {
relativeParetoWeightUpperBound <- uniroot(
function(x) {calculateDiffLCRiskSharingAndAutarky(
x,
relativeParetoWeightsGridPoints,
delta,
sigma,
aggregateIncome,
householdValueLCRiskSharing,
villageValueLCRiskSharing,
incomeTransitionProbVec,
householdAutarkyValue,
villageAutarkyValue,
numHouseholds
)$villageDiffLCRiskSharingAndAutarky},
c(min(relativeParetoWeightsGridPoints), max(relativeParetoWeightsGridPoints)),
tol = 1e-10,
maxiter = 300
)$root
} else {
relativeParetoWeightUpperBound <- max(relativeParetoWeightsGridPoints)
}
relativeParetoWeightBounds[incomeStateIndex, 1] <- relativeParetoWeightLowerBound
relativeParetoWeightBounds[incomeStateIndex, 2] <- relativeParetoWeightUpperBound
}
if (iteration == iterationLimit) {
print("Reached the maximum limit of iterations!")
}
return(list(
householdValueLCRiskSharing = householdValueLCRiskSharing,
villageValueLCRiskSharing = villageValueLCRiskSharing,
relativeParetoWeightBounds = relativeParetoWeightBounds))
}solveLCRiskSharing <- function(
delta,
sigma,
punishment,
householdIncomeGridPoints,
villageIncomeGridPoints,
incomeTransitionMatrix,
incomeGridPointsMatrix,
numIncomeStates,
numHouseholds,
numRelativeParetoWeights = 2000,
iterationLimit = 100,
diffLimit = 1e-8
) {
aggregateIncomeGridPoints <- (
incomeGridPointsMatrix[, 1] + incomeGridPointsMatrix[, 2] * (numHouseholds - 1)
)
autarkyValueMatrix <- cbind(
calculateAutarkyValue(
householdIncomeGridPoints,
sigma,
delta,
punishment,
incomeTransitionMatrix
),
calculateAutarkyValue(
villageIncomeGridPoints,
sigma,
delta,
punishment,
incomeTransitionMatrix
)
)
relativeParetoWeightsGridPoints <- getRelativeParetoWeightsGridPoints(
sigma,
punishment,
householdIncomeGridPoints,
villageIncomeGridPoints,
numRelativeParetoWeights
)
consumptionOnRelativeParetoWeightGrid <- matrix(
NA, nrow = numIncomeStates, ncol = numRelativeParetoWeights
)
for (incomeStateIndex in seq_along(aggregateIncomeGridPoints)) {
for (relativeParetoWeightIndex in seq_along(relativeParetoWeightsGridPoints)) {
consumptionOnRelativeParetoWeightGrid[
incomeStateIndex,
relativeParetoWeightIndex
] <- calculateHouseholdConsumption(
aggregateIncomeGridPoints[incomeStateIndex],
relativeParetoWeightsGridPoints[relativeParetoWeightIndex],
numHouseholds,
sigma
)
}
}
valueFullRiskSharing <- calculateValueFullRiskSharing(
incomeTransitionMatrix,
aggregateIncomeGridPoints,
delta,
sigma,
autarkyValueMatrix,
consumptionOnRelativeParetoWeightGrid,
numRelativeParetoWeights,
numHouseholds
)
valueLCRiskSharing <- calculateValueLCRiskSharing(
valueFullRiskSharing,
consumptionOnRelativeParetoWeightGrid,
aggregateIncomeGridPoints,
incomeTransitionMatrix,
autarkyValueMatrix,
relativeParetoWeightsGridPoints,
numRelativeParetoWeights,
delta,
sigma,
numIncomeStates,
numHouseholds,
iterationLimit,
diffLimit
)
return(valueLCRiskSharing)
}sigma <- 1.0
punishment <- 0.0
delta <- 0.8
sigma <- 1
numHouseholds <- 2
numIncomeStates <- 3
incomeGridPointsHH1 <- c(0.353, 0.5, 0.647)
incomeGridPointsHH2 <- 1 - incomeGridPointsHH1
incomeGridPointsMatrix <- cbind(incomeGridPointsHH1, incomeGridPointsHH2)
aggregateIncomeGridPoints <- incomeGridPointsHH1 + incomeGridPointsHH2
incomeTransitionMatrix <- matrix(1 / 3, nrow = 3, ncol = 3)LCRiskSharingNoStorageResult <- solveLCRiskSharing(
delta,
sigma,
punishment,
incomeGridPointsHH1,
incomeGridPointsHH2,
incomeTransitionMatrix,
incomeGridPointsMatrix,
numIncomeStates,
numHouseholds,
numRelativeParetoWeights = 10000,
iterationLimit = 1000,
diffLimit = 1e-8
)set.seed(123)
simSize <- 5000
incomeSeq <- sample(seq(1, numIncomeStates), size = simSize, replace = TRUE)
relativeParetoWeightSeqNoStorage <- c(1)
consHH1SeqNoStorage <- c()
for (i in seq(1, simSize)) {
relativeParetoWeightSeqNoStorage <- c(
relativeParetoWeightSeqNoStorage,
relativeParetoWeightSeqNoStorage[i] %>%
pmax(
LCRiskSharingNoStorageResult$relativeParetoWeightBounds[
incomeSeq[i],
1
]
) %>%
pmin(
LCRiskSharingNoStorageResult$relativeParetoWeightBounds[
incomeSeq[i],
2
]
)
)
consHH1SeqNoStorage <- c(
consHH1SeqNoStorage,
calculateHH1Consumption(
(aggregateIncomeGridPoints[incomeSeq[i]]),
relativeParetoWeightSeqNoStorage[i + 1],
sigma
)
)
}
simulationResultNoStorage <- tibble(
period = seq(1, simSize),
income = incomeSeq,
consHH1 = consHH1SeqNoStorage
) %>%
mutate(
previousIncome = lag(income),
changeIncome = paste(as.character(previousIncome), '->', as.character(income))
)4.2 Simulation result
The figure shows the consumption evolution over time. Unlike the model with storage, the model without storage results in a consumption pattern which takes only a finite number of values in the steady state.
simulationResultNoStorage %>%
filter(period >= 100) %>%
ggplot() +
geom_point(aes(x = period, y = consHH1), size = 0.3) +
xlab('Time period') +
ylab('HH1 consumption') +
theme_classic()
4.3 Comparison with the case with storage
The table below compares the consumption means, standard deviations, and the average utilities over time. Due to the storage with positive return on saving, the mean consumption is higher in the model with storage. However, due to the existence of storage, the autarky values are higher, and hence the risk-sharing functionality is degraded. This results in more variable consumption (second row) and hence, lower welfare (third row) in the model with storage. Notice that, this kind of reduced welfare is achieved only for medium-level returns on saving. If the return is too low, saving is not used and there is no welfare impact, and if the return is sufficiently high, the benefit from intertemporal consumption smoothing dominates and hence there is a welfare improvement.
simulationResult <- readRDS(
file.path('IntermediateData/simulationResult.rds')
)summaryTable <- cbind(
c(
(simulationResult %>% filter(period >= 100) %>% .$consHH1) %>% mean,
(simulationResult %>% filter(period >= 100) %>% .$consHH1) %>% sd,
log(simulationResult %>% filter(period >= 100) %>% .$consHH1) %>% mean
),
c(
(simulationResultNoStorage %>% filter(period >= 100) %>% .$consHH1) %>% mean,
(simulationResultNoStorage %>% filter(period >= 100) %>% .$consHH1) %>% sd,
log(simulationResultNoStorage %>% filter(period >= 100) %>% .$consHH1) %>% mean
)
)
rownames(summaryTable) <- c("Mean Cons.", "SD Cons", "Mean Util.")
colnames(summaryTable) <- c("With storage", "Without storage")
summaryTable %>%
kbl(digits = 4) %>%
kable_classic()| With storage | Without storage | |
|---|---|---|
| Mean Cons. | 0.5009 | 0.4992 |
| SD Cons | 0.0863 | 0.0717 |
| Mean Util. | -0.7067 | -0.7053 |