We first clear the workspace using rm(list = ls())and
then include all packages we need. If a package is missing in your R
distribution (which is quite likely initially), just use
install.packages("package_name") with the respective
package name to install it on your system. If you execute the code in
the file install_packages.R, then all necessary packages
will be installed into your R distribution. If the variable
export_graphs is set to TRUE, then the graphs
will be exported as pdf-files. In addition, we define a set of colors
here to make graphs look more beautiful.
rm(list = ls())
library(reshape2)
library(base)
library(ggplot2)
library(grid)
library(scales)
library(stringr)
library(tidyverse)
# should graphs be exported to pdf
export_pdf <- FALSE
# define some colors
mygreen <- "#00BA38"
myblue <- "#619CFF"
myred <- "#F8766D"
Before being able to simulate the model, we have to set parameters. To do so, we use five data moments. The share of labor compensation in GDP, the share of gross capital formation (investment) in GDP, the gross return to capital and the growth rates of total labor input and labor productivity. From these, we can pin down the parameters of the Ramsey model for the US economy. Having set the parameters, we can determine the balanced growth path, i.e. the steady state of the capital intensity.
# data moments
lab_share = 0.381
inv_share = 0.245
gross_ret = 0.114
# parameter choice
n <- 0.0122
h <- 0.0169
g <- (1+n)*(1+h)-1
alpha <- lab_share
delta <- inv_share*gross_ret/alpha - g
beta <- (1+g)/(1 + gross_ret - delta)
gamma <- 0.5
# steady state capital intensity under old beta
k_old <- (alpha/((1+g)/beta - 1 + delta))^(1/(1-alpha))
c_old <- k_old^alpha - (g+delta)*k_old
In the lecture, we discussed a guessing algorithm for solving for the saddle path of the Ramsey economy. This algorithm is implemented in function `solve_ramsey’.
The function takes as input the relevant parameters
g, alpha, delta, gamma, beta of the Ramsey model. In
addition, we have to specify the level of the capital intensity
k0 at which we want to start the transition path towards
the balanced growth path.
Within the function, we have to specify a set of numerical
parameters. T is the maximum number of model periods we
will simulate in the algorithm. itermax is the maximum
number of guesses for \(c_0\) we will
use to solve for the saddle path. tol is the level of
tolerance. If the distance of \(k_t\)
and \(c_t\) to its steady state values
is smaller than tol, we say that the economy has
converged.
Next we determine the steady state capital intensity of the Ramsey and the corresponding consumption level. After that, we can determine the relevant bounds for \(c_0\) in the guessing algorithm. If \(k_0 < k^*\), then \(c_0\) must lie in the interval \([0, c^*]\). Otherwise, a good guess for the relevant interval is \([c^*, k_0^{\alpha} + (1-delta)*k_0]\).
With this setup, we can start the guessing algorithm discussed in class. We set the guess for initial consumption \(c_0 = \frac{c_{min} + c_{max}}{2}\) and simulate the economy forward until:
The capital-consumption path has converged to the steady state of the capital intensity. We define convergence as \[ \varepsilon = \sqrt{\left(\frac{c_t}{c^*}-1\right)^2 + \left(\frac{k_t}{k^*}-1\right)^2} < tol. \]
The capital intensity exceeds \(k^*\), which means that initial consumption was too low. We therefore set \(c_{min}\) to the last initial guess for consumption.
Consumption exceeds \(c^*\), which means that initial consumption was too high. We therefore set \(c_{max}\) to the last initial guess for consumption.
Note that the above algorithm describes a situation with \(k_0 \leq k^*\). In case of \(k_0 > k^*\), the economy approaches the balanced growth path from above. Hence, situations 2. and 3. are flipped.
Our function returns three arguments: the best guess for initial consumption \(c_0\), the convergence level \(\varepsilon\) as well as the number of iteration \(i\) needed to find the initial consumption level.
solve_ramsey <- function(k0, g, alpha, delta, gamma, beta) {
# set numerical parameters
T <- 10000 # number of simulation periods
itermax <- 1000 # maximum number of iterations
tol <- 1e-6 # tolerance level for difference to steady state
# determine steady state capital intensity and consumption
kstar <- (alpha/((1+g)/beta - 1 + delta))^(1/(1-alpha))
cstar <- kstar**alpha - (g+delta)*kstar
# set starting interval depending on k0
if(k0 <= kstar) {
c_min <- 0
c_max <- cstar
} else {
c_min <- cstar
c_max <- k0^alpha + (1-delta)*k0
}
# start iteration process
for (i in 1:itermax) {
# start a capital and consumption path
k <- k0
c <- (c_min+c_max)/2
# now simulate the economy forward
for (t in 1:T) {
# calculate difference to steady state
epsilon <- sqrt((c[t]/cstar-1)^2 + (k[t]/kstar-1)^2)
# if difference is small enough, that's it, we return
if(epsilon < tol) {
return(c(c[1], epsilon, i))
}
# determine next period's capital stock and consumption
k[t+1] <- max((k[t]**alpha - c[t] + (1-delta)*k[t])/(1+g), 1e-4)
c[t+1] <- (beta*(1 + alpha*k[t+1]^(alpha-1) - delta)/(1+g))^gamma*c[t]
# decision rule for capital
if(k0 <= kstar & k[t+1] > kstar | k0 > kstar & c[t+1] < cstar) {
c_min <- c[1]; break;
}
# if consumption is too large, then initial consumption level was too small
if(k0 <= kstar & c[t+1] > cstar | k0 > kstar & k[t+1] < kstar) {
c_max <- c[1]; break;
}
}
}
# in any case, if you end up here, return c[1] and epsilon
return(c(c[1], epsilon, i))
}
Once we have determined the initial consumption level \(c_0\) on the saddle path, we can use
standard forward simulation to determine the entire path to the steady
state. This is done in the function ramsey, which follows
the same logic as the simulation of the
Solow growth model. We again simulate the economy from some date
\(T_0 < 0\) to some date \(T_1 > 0\). We assume that the economy is
in steady state with capital intensity k0 until date \(t = -1\). The change in \(\beta\) is assumed to happen at date \(t = 0\). The economy then starts on a
transition path with initial capital intensity k0. We again
use a function ind that transforms dates into array
indices, since array indices start at \(1\).
For forward simulation of the economy, we use the two equations \[ k_{t+1} = \max\left[\frac{k_t^\alpha - c_t + (1-\delta)k_t}{1+g}, 10^{-4}\right] \] and \[ c_{t+1} = \left[\beta\cdot\frac{1 + \alpha k_{t+1}^{\alpha-1} - \delta}{1+g}\right]^\gamma\cdot c_t. \] Note that we restrict \(k_{t+1}\) to not fall below a value of \(10^{-4}\), such that we don’t run into calculation errors on the marginal product of capital. Furthermore, if the economy is reasonably close to the steady state – according to the same criterion we used when guessing the initial point of the transition path – then we map it onto the steady state capital and consumption level. In doing so, we avoid that – owing to numerical inaccuracies – the path for capital diverges.
ramsey <- function(T0, T1, c0, k0, g, alpha, delta, gamma, beta) {
# set numerical parameters
tol <- 1e-6
# determine steady state capital intensity and consumption
kstar <- (alpha/((1+g)/beta - 1 + delta))^(1/(1-alpha))
cstar <- kstar**alpha - (g+delta)*kstar
# start a capital and consumption path
k <- k0
c <- 0
# assume economy was in steady state prior to date 0
k[ind(T0):ind(0)] <- k0
c[ind(T0):ind(0)] <- k0^alpha - (g+delta)*k0
# then start to simulate new economy path
c[ind(0)] <- c0
# now simulate the economy forward
for (t in ind(0):ind(T1-1)) {
# determine next period's capital stock and consumption
k[t+1] <- max((k[t]**alpha - c[t] + (1-delta)*k[t])/(1+g), 1e-4)
c[t+1] <- (beta*(1 + alpha*k[t+1]^(alpha-1) - delta)/(1+g))^gamma*c[t]
# move to steady state level when you are almost there
epsilon <- sqrt((c[t+1]/cstar-1)^2 + (k[t+1]/kstar-1)^2)
if(epsilon < tol) {
k[(t+1):ind(T1)] = kstar
c[(t+1):ind(T1)] = cstar
break
}
}
# calculate other macro statistics
y <- k^alpha
ir <- (g + delta)*k
i <- y - c
s <- i/y
# return a data frame with macro path
res <- data.frame(year=c(T0:T1), k, y, c, i, ir, s)
return(res)
}
# indicator management function
ind <- function(t) {
return(t + 1 + abs(T0))
}
When we want to simulate an increase in the discount factor \(\beta\), we first have to calculate the balanced growth path capital intensity \(k_{old}\) under the original discount factor. Then, we can set the parameters for the simulation process. We want to increase the discount factor to \(\beta = 0.99\), simulate 50 periods prior to the change in \(\beta\) and then let the economy run for 200 periods starting from \(t=0\) onward.
With all parameters at hand, we can solve for the initial consumption
level \(c_0\) that defined the
beginning of the transition path to the new balanced growth path with a
higher discount factor, i.e. we trace out the beginning of the saddle
path. The function solve_ramsey returns \(c_0\) as a first argument of the vector
res. The second argument is the accuracy level with which
we hit the balanced growth path and the third the number of iterations
we needed to calculate \(c_0\).
Once \(c_0\) is known, we can
simulate the entire transition path with the function
ramsey. This function takes the parameters for the
simulation periods, the initial capital level \(k_{old}\) and the initial level of
consumption stored in res[1]. It returns a data frame that
contains the entire path of the economy. We store this data frame in
transition and can use it to create the relevant plots.
# steady state capital intensity under old beta
k_old <- (alpha/((1+g)/beta - 1 + delta))^(1/(1-alpha))
c_old <- k_old^alpha - (g+delta)*k_old
# set new value for beta and simulation periods
beta <- 0.99
T0 <- -50
T1 <- 200
# solve transition path
res <- solve_ramsey(k_old, g, alpha, delta, gamma, beta)
# simulate the model for 250 periods
transition <- ramsey(T0, T1, res[1], k_old, g, alpha, delta, gamma, beta)
We first plot the dynamics of the capital intensity. When the
discount factor increases, so does the long-run capital intensity.
Hence, the economy is put on a transition path, where the capital stock
increase successively and converges toward the new long-run capital
intensity. This looks pretty much like our analysis of the
Solow growth model.
myplot <- ggplot(data = transition) +
geom_hline(yintercept=transition$k[ind(T0)], color=myred, linetype="dashed", linewidth=1) +
geom_hline(yintercept=transition$k[ind(T1)], color="#00BA38", linetype="dashed", linewidth=1) +
geom_line(aes(x=year, y=k), color="darkblue", linewidth=1) +
coord_cartesian(xlim=c(T0, T1), ylim=c(6, 14)) +
scale_x_continuous(breaks=seq(T0, T1, 25), expand=c(0, 0)) +
labs(x = "Year t",
y = "Capital Intensity") +
theme_bw()
# print the plot
print(myplot)
What is different compared to the Solow model, however, is the
dynamics of the savings rate. Recall that the savings rate is an
exogenous object in the Solow model, and in our analysis it just jumped
from one date to the other. In the Ramsey model, the savings rate is an
endogenous object that is determined by intertemporal optimization. We
can plot the dynamics of the savings rate in our model, which we can
calculate from relating investment to GDP. We can clearly see that the
savings rate is a time varying object. It is initially higher than in
the balanced growth path and then converges to the balanced growth
savings rate (which is the same as we would see in the Solow model).
myplot <- ggplot(data = transition) +
geom_hline(yintercept=(g+delta)/transition$k[ind(T0)]^(alpha-1), color=myred, linetype="dashed", linewidth=1) +
geom_hline(yintercept=(g+delta)/transition$k[ind(T1)]^(alpha-1), color="#00BA38", linetype="dashed", linewidth=1) +
geom_line(aes(x=year, y=s), color="darkblue", linewidth=1) +
coord_cartesian(xlim=c(T0, T1), ylim=c(0, 0.5)) +
scale_x_continuous(breaks=seq(T0, T1, 25), expand=c(0, 0)) +
labs(x = "Year t",
y = "Savings Rate s") +
theme_bw()
# print the plot
print(myplot)
Finally, we can again look at the dynamics of GDP and its
components. We get a very similar picture as in our analysis of the
Solow growth model. Household initially cut back on their
consumption to invest above the replacement investment rate. This slowly
increases the capital stock. In the long run, investment again equals
replacement investment, but at a higher capital intensity.
myplot <- ggplot(data = transition) +
geom_hline(yintercept = c_old, color="#00BA38", linetype="dashed", linewidth=0.5) +
geom_ribbon(aes(x=year, ymin=0, ymax=c, fill= "1c", color="1c") , alpha=0.4) +
geom_ribbon(aes(x=year, ymin=c, ymax=c+i-ir, fill= "3di", color="3di"), alpha=0.4) +
geom_ribbon(aes(x=year, ymin=c+i-ir, ymax=y, fill= "2ir", color="2ir") , alpha=0.4) +
geom_line(aes(x=year, y=y), color="darkblue", linewidth=1) +
coord_cartesian(xlim=c(T0, T1), ylim=c(0, 3)) +
scale_x_continuous(breaks=seq(T0, T1, 25), expand=c(0, 0)) +
labs(x = "Year t",
y = "GDP and its components per\n effective unit of labor") +
scale_fill_manual(breaks = c("1c", "2ir", "3di"), name = "",
labels = c("Consumption", "Replacement Investment", "Capital Augmenting Inv."),
values = c(mygreen, myblue, myred)) +
scale_color_manual(breaks = c("1c", "2ir", "3di"),
values = c(mygreen, myblue, myred)) +
guides(colour = "none") +
theme_bw() +
theme(legend.position="bottom")
# print the plot
print(myplot)
Now we want to perform the same exercise as above, but with a value of \(\gamma = 0.1\). We therefore use the same code as above.
# set new gamma and old beta
gamma <- 0.1
beta <- (1+g)/(1 + gross_ret - delta)
# steady state capital intensity under old beta
k_old <- (alpha/((1+g)/beta - 1 + delta))^(1/(1-alpha))
c_old <- k_old^alpha - (g+delta)*k_old
# set new value for beta and simulation periods
beta <- 0.99
T0 <- -50
T1 <- 200
# solve transition path
res <- solve_ramsey(k_old, g, alpha, delta, gamma, beta)
# simulate the model for 250 periods
transition <- ramsey(T0, T1, res[1], k_old, g, alpha, delta, gamma, beta)
Compared to the case of \(\gamma =
0.5\), the dynamics of the capital stock is much slower.
# Plot dynamics of the capital stock
myplot <- ggplot(data = transition) +
geom_hline(yintercept=transition$k[ind(T0)], color=myred, linetype="dashed", linewidth=1) +
geom_hline(yintercept=transition$k[ind(T1)], color="#00BA38", linetype="dashed", linewidth=1) +
geom_line(aes(x=year, y=k), color="darkblue", linewidth=1) +
coord_cartesian(xlim=c(T0, T1), ylim=c(6, 14)) +
scale_x_continuous(breaks=seq(T0, T1, 25), expand=c(0, 0)) +
labs(x = "Year t",
y = "Capital Intensity") +
theme_bw()
# print the plot
print(myplot)
The reason for this can be found by looking at the dynamics of
the savings rate. While under \(\gamma =
0.5\), the savings rate jumped above its long-run level, because
households wanted to accumulate savings quickly to match their new
discount factor, the savings rate is now much more slow moving over
time. It starts at a value below the long-run level and then slowly
increases over time. This keeps the consumption path much smoother than
under \(\gamma = 0.5\).
# Plot dynamics of the savings rate
myplot <- ggplot(data = transition) +
geom_hline(yintercept=(g+delta)/transition$k[ind(T0)]^(alpha-1), color=myred, linetype="dashed", linewidth=1) +
geom_hline(yintercept=(g+delta)/transition$k[ind(T1)]^(alpha-1), color="#00BA38", linetype="dashed", linewidth=1) +
geom_line(aes(x=year, y=s), color="darkblue", linewidth=1) +
coord_cartesian(xlim=c(T0, T1), ylim=c(0, 0.5)) +
scale_x_continuous(breaks=seq(T0, T1, 25), expand=c(0, 0)) +
labs(x = "Year t",
y = "Savings Rate s") +
theme_bw()
# print the plot
print(myplot)
We can see the same picture when plotting GDP and its
components. The initial drop in consumption is much smaller than under
\(\gamma = 0.5\). This also means that
capital augmenting investment is much smaller initially and therefore
GDP rises more slowly over time.
# Plot GDP and its components
myplot <- ggplot(data = transition) +
geom_hline(yintercept = c_old, color="#00BA38", linetype="dashed", linewidth=0.5) +
geom_ribbon(aes(x=year, ymin=0, ymax=c, fill= "1c", color="1c") , alpha=0.4) +
geom_ribbon(aes(x=year, ymin=c, ymax=c+i-ir, fill= "3di", color="3di"), alpha=0.4) +
geom_ribbon(aes(x=year, ymin=c+i-ir, ymax=y, fill= "2ir", color="2ir") , alpha=0.4) +
geom_line(aes(x=year, y=y), color="darkblue", linewidth=1) +
coord_cartesian(xlim=c(T0, T1), ylim=c(0, 3)) +
scale_x_continuous(breaks=seq(T0, T1, 25), expand=c(0, 0)) +
labs(x = "Year t",
y = "GDP and its components per\n effective unit of labor") +
scale_fill_manual(breaks = c("1c", "2ir", "3di"), name = "",
labels = c("Consumption", "Replacement Investment", "Capital Augmenting Inv."),
values = c(mygreen, myblue, myred)) +
scale_color_manual(breaks = c("1c", "2ir", "3di"),
values = c(mygreen, myblue, myred)) +
guides(colour = "none") +
theme_bw() +
theme(legend.position="bottom")
# print the plot
print(myplot)