```python
from scipy import optimize
import math as math
import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
import ipywidgets as wd
```

# A simple CGE model
In the first section of this code we construct a simple Computable General Equilibruim (CGE) model and calibrate it to an observed equilibrium. In the second section of this code we examine the properties of the model and the effect of certain shocks.

## Building a CGE model - a crash course in CGE-modelling
The building blocks of a CGE model are:
1. Model equations describing an economy
2. Behavioural parameters
3. Variable values representing an initial equilibrium

Below we will construct a simple CGE model in three steps:
1. We will define the model's equations. The model is a simple open economy with producers, consumers and a world market
2. We will discuss the behavioural parameters
3. We will calibrate the models unknown behavioural parameters such that the model can replicate the initial equilibrium

### 1. Model equations
Our model is a slightly modified version of a model from the course in applied CGE-modelling. The model consists of producers who buy labour and produce a single domestic consumption good, cosumers who sell labour to domestic producers and buy both domestic and foregin consumption goods, and a world market where producers can sell sell domestic goods and consumers can buy foreign goods.
#### Domestic producers
The domestic producers are in full competition with each other and all have the CES-production function:

$$ Y = \left[ \mu^\frac{1}{E_Y} L^\frac{E_Y-1}{E_Y} \right]^\frac{E_Y}{E_Y-1} $$

They maximize their profits:

$$ \pi = p_dY - wL $$

This implies the following labour demand and zero profit condition, i.e. the price of the domestic consumption good is equal to the avg. costs of producing it.

$$ 
\begin{align}
L = \mu  \left( \frac{w}{p_d} \right)^{-E_Y} Y \tag{1}
\end{align}
$$

$$ 
\begin{align}
p_d = \frac{w L}{Y} \tag{2} 
\end{align}
$$

Below we define the two equations:


```python
## Labour demand
def E1_L_D(mu,w,p_d,E_Y,Y) :
    """ Funtion defining the producers' demand for labour. """
    return mu * (w/p_d)**(-E_Y) * Y

## Zero profit condition - expressed as a function of p_d
def E2_p_d(Y,w,L) :
    """  Function defining the zero profit condidtion. """
    return w*L / Y
```

#### Domestic consumers
The model has N domestic consumers who have a CES utility function and derive utility from the consumption of domestic and foreign goods:

$$ U(C_d,C_f)=\left[ \gamma_d^\frac{1}{E_C} C_d^\frac{E_C-1}{E_C} + \gamma_f^\frac{1}{E_C} C_f^\frac{E_C-1}{E_C} \right]^\frac{E_C}{E_C-1} $$

Their budget constaint is:

$$
\begin{align}
p_d C_d + p_f C_f = w N \tag{3}
\end{align}
$$

The consumers maximize their utility subject to their budget constraint. This implies the following demand functions:

$$
\begin{align}
C_d = \gamma_d \left(\frac{p_d}{P_C}\right)^{-E_C} \frac{wN}{P_C} \tag{4}
\end{align}
$$

$$
\begin{align}
C_f = \gamma_d \left(\frac{p_f}{P_C}\right)^{-E_C} \frac{wN}{P_C} \tag{5}
\end{align}
$$

$P_C$ is a CES-price index.


```python
## The budget of the households constraint
def E3_wN(p_d,C_d,p_f,C_f) :
    """ Budget constraint of the households. """
    return p_d*C_d+p_f*C_f

## The demand for the domestic consumption good
def E4_C_d(gamma_d,p_d,P_C,E_C,w,N) :
    """ Demand for domestic consumption."""
    return gamma_d * (p_d/P_C)**(-E_C) *w*N/P_C

## The demand for the foreign consumption good
def E5_C_f(gamma_f,p_f,P_C,E_C,w,N) :
    """ Demand for foregin consumption good."""
    return gamma_f * (p_f/P_C)**(-E_C) *w*N/P_C
```

#### The labour market
We assume that everybody works. Thus the labour market becomes:

$$
\begin{align}
L = N \tag{6}
\end{align}
$$


```python
## The labour market
def E6_L_S(N) :
    """ Labour market. """
    return N
```

#### The goods market
In equilibria the entire domestic production is either consumed by domestic consumers or exported:

$$ 
\begin{align}
Y = C_d + X \tag{7}
\end{align}$$


```python
## The goods market
def E7_Y(C_d,X) :
    """ The goods market. """
    return C_d + X
```

#### Foregin trade
We model the foreign trade using an Armington's approach, where the export depends on the relative price of foreign and domestic goods:

$$ 
\begin{align}
X = \phi \left(\frac{p_d}{p_f}\right)^{-E_X} \tag{8}
\end{align}
$$

We use the foreign price as a numéraire:

$$
\begin{align}
p_f = 1 \tag{9}
\end{align}
$$


```python
## Armington's approach to foreign trade
def E8_X(phi,p_d,p_f,E_X) :
    """ Foreign trade """
    return phi*(p_d/p_f)**(-E_X)

## Setting the foreign price as a numéraire
def E9_p_f() :
    """ The price of the foreign good """
    return 1
```

#### The model
Equation (1)-(9) represents the nine equations in our model. We treat the nine variables $Y$, $L$, $w$, $C_d$, $C_f$, $p_d$, $p_f$, $P_C$, and $X$ as endogenous.

The remaining parameters are $N$, $\mu$, $\gamma_d$, $\gamma_f$, $\phi$, $E_Y$, $E_C$ and $E_X$.

### 2. Behavioural parameters
The model has 8 unknown parameters: $N$, $\mu$, $\gamma_d$, $\gamma_f$, $\phi$, $E_Y$, $E_C$ and $E_X$. $\mu$, $\gamma_d$, $\gamma_f$, $\phi$, $E_Y$, $E_C$ and $E_X$ are the so-called behavioural parameters in the CES-functions and in Armington's approach to foreign trade. In the original assignment from the course in applied CGE-modelling, it was given that $E_Y \equiv 2.0$, $E_C \equiv 0.5$, and $E_C \equiv 5.0$.

The remaining behavioural parameters and $N$ are determined by calibrating the model to the initial equilibrium, e.g. setting $\phi$ in equation $(8)$ such that the equation is consistent with the initial values of $X$, $p_d$, and $p_f$.


```python
E_Y = 2.0
E_C = 0.5
E_X = 5.0
```

### 3. Initial dataset and calibration
To model is calibrated to an initial equilibrium described by the IO-table:

| I/O | $PS$   | $PC$  | $X$   |
| --- | --:  | --: | --: |
| $PS$  | 0    | 800 | 200 |
| $M$   | 0    | 200 |   0 |
| $w$ | 1000 | 0   |   0 |

PS is the private sector, $M$ is imports, $w$ is wages, PC is private consumption, and $X$ is exports. The rows describe input and the columns output. The private sector e.g. uses 1.000 units of wage, and outputs 800 units of goods used for domestic consumption. The table represents a simple National Account in current prices.


```python
rows    = ['PS','M' ,'w']
columns = ['PS','PC','X']

data = [(0    , 800 , 200),
        (0    , 200 ,   0) ,
        (1000 , 0   ,   0) ]

IO = pd.DataFrame(data,columns=columns,index=rows)
```

#### Initialising the variables
All prices are simply set to 1 and interpretated as price indices. The ammount of labour, the total prodcution etc. are found using the IO-table. The IO-table is in current prices so we devide by the price indices.


```python
## Endogenous variables in initial equilibrium
# Prices
p_d0 = 1
p_f0 = 1
w0   = 1
P_C0 = 1

# Amounts
L0   = IO.loc['w']['PS']/w0
Y0   = IO.loc['w']['PS']/p_d0
X0   = IO.loc['PS']['X']/p_d0
C_d0 = IO.loc['PS']['PC']/p_d0
C_f0 = IO.loc['M']['PC']/p_f0

## List of initial values
ini_list = [Y0,L0,w0,C_d0,C_f0,p_d0,p_f0,P_C0,X0]
```

We also save the initial equilibrium in a dataframe.


```python
## Saving initial equilibrium in a dictionary
ini_eq = [{"Shock" : "Initial eq.", 
            "$L$"   : "%1.f"  % L0,
            "$Y$"   : "%0.2f" % Y0,
            "$w$"   : "%0.3f" % w0,
            "$C_d$" : "%0.2f" % C_d0,
            "$p_d$" : "%0.3f" % p_d0,
            "$C_f$" : "%0.2f" % C_f0,
            "$p_f$" : "%0.3f" % p_f0,
            "$P_C$" : "%0.3f" % P_C0,
            "$X$"   : "%0.2f" % X0}]

## Converting the dictionary to a data-frame
ini_eq = pd.DataFrame(ini_eq)

## Rearranging the columns of the dataframe
cols = ['Shock','$L$','$Y$','$X$','$C_d$','$C_f$','$w$','$p_d$','$p_f$','$P_C$']
ini_eq = ini_eq[cols]

## Printing the inital equilibrium
ini_eq
```




<div>
<style scoped>
    .dataframe tbody tr th:only-of-type {
        vertical-align: middle;
    }

    .dataframe tbody tr th {
        vertical-align: top;
    }

    .dataframe thead th {
        text-align: right;
    }
</style>
<table border="1" class="dataframe">
  <thead>
    <tr style="text-align: right;">
      <th></th>
      <th>Shock</th>
      <th>$L$</th>
      <th>$Y$</th>
      <th>$X$</th>
      <th>$C_d$</th>
      <th>$C_f$</th>
      <th>$w$</th>
      <th>$p_d$</th>
      <th>$p_f$</th>
      <th>$P_C$</th>
    </tr>
  </thead>
  <tbody>
    <tr>
      <th>0</th>
      <td>Initial eq.</td>
      <td>1000</td>
      <td>1000.00</td>
      <td>200.00</td>
      <td>800.00</td>
      <td>200.00</td>
      <td>1.000</td>
      <td>1.000</td>
      <td>1.000</td>
      <td>1.000</td>
    </tr>
  </tbody>
</table>
</div>



#### Calibrating
The next step is to calibrate $N$, $\mu$, $\gamma_d$, $\gamma_f$, $\phi$, i.e. setting their values such that the model replicates the initial equilibrium. In the case of e.g. $\mu$ this is done by solving the labour demand function, equation (1), for $\mu$ with the inital values inserted:

$$ 
\begin{align}
L_0 = \mu  \left( \frac{w_0}{p_{d0}} \right)^{-E_Y} Y_0
\end{align}
$$

We do this numerically below.


```python
## Calibrating population size
def f_L_S(N) :
    return E6_L_S(N) - L0

solution = optimize.root(f_L_S, (0) )
N = np.asscalar(solution.x)

## Calibrating mu
def f_mu(mu) :
    return E1_L_D(mu,w0,p_d0,E_Y,Y0) - L0

solution = optimize.root(f_mu, (0) )
mu = np.asscalar(solution.x)

## Calibrating gamma_d and gamma_f
def f_gamma_d(gamma_d) :
    return E4_C_d(gamma_d,p_d0,P_C0,E_C,w0,N) -C_d0

solution = optimize.root(f_gamma_d, (0) )
gamma_d = np.asscalar(solution.x)

def f_gamma_f(gamma_f) :
    return E5_C_f(gamma_f,p_f0,P_C0,E_C,w0,N) -C_f0

solution = optimize.root(f_gamma_f, (0) )
gamma_f = np.asscalar(solution.x)

## Calibrating phi
def f_phi(phi) :
    return E8_X(phi,p_d0,p_f0,E_X) - X0

solution = optimize.root(f_phi, (0) )
phi = np.asscalar(solution.x)

## Printing the results
print("N =",N,", mu =",mu,",","gamma_d =",gamma_d,",","gamma_f =",gamma_f,",","phi =",phi)
```

    N = 1000.0 , mu = 1.0 , gamma_d = 0.8 , gamma_f = 0.2 , phi = 200.0


#### Solving the model
The last step of building our model is to solve it and check if it replicates the initial equilibrium. If it does then our model works. This is called a zero-shock in CGE-modelling.

The model is solved as a system of nine equations, equation $(1)$-$(9)$, and nine unknowns, $Y$, $L$, $w$, $C_d$, $C_f$, $p_d$, $p_f$, $P_C$, and $X$:

$$ 
\begin{align}
L = \mu  \left( \frac{w}{p_d} \right)^{-E_Y} Y \tag{1}
\end{align}
$$

$$ 
\begin{align}
p_d = \frac{w L}{Y} \tag{2} 
\end{align}
$$

$$
\begin{align}
p_d C_d + p_f C_f = w N \tag{3}
\end{align}
$$

$$
\begin{align}
C_d = \gamma_d \left(\frac{p_d}{P_C}\right)^{-E_C} \frac{wN}{P_C} \tag{4}
\end{align}
$$

$$
\begin{align}
C_f = \gamma_d \left(\frac{p_f}{P_C}\right)^{-E_C} \frac{wN}{P_C} \tag{5}
\end{align}
$$

$$
\begin{align}
L = N \tag{6}
\end{align}
$$

$$ 
\begin{align}
Y = C_d + X \tag{7}
\end{align}
$$

$$ 
\begin{align}
X = \phi \left(\frac{p_d}{p_f}\right)^{-E_X} \tag{8}
\end{align}
$$

$$
\begin{align}
p_f = 1 \tag{9}
\end{align}
$$

We solve the model by defining a function called CGEsolve. The function defines the system of equations and solves it. This functionh will become handy later, when we examine the properties of the model.


```python
## Function solving the CGE model taking parameters as inputs
def CGEsolve(N,mu,gamma_d,gamma_f,phi,E_Y,E_C,E_X,ini_list,status='yes') :
    """
    This function defines and solves the CGE-model as a function of its parameters and
    a list of inital values for the solver.
    
    status = 'yes' prints a summary of the solver results
    """
    ## Defining the CGE model as a system of nine equations and nine unknowns
    def CGEmodel(variables) :
        ## Defining variables
        (Y,L,w,C_d,C_f,p_d,p_f,P_C,X) = variables

        ## Defining equations
        EQ_L_D = E1_L_D(mu,w,p_d,E_Y,Y) - L
        EQ_p_d = E2_p_d(Y,w,L) - p_d
        EQ_wN  = E3_wN(p_d,C_d,p_f,C_f) - w*N
        EQ_C_d = E4_C_d(gamma_d,p_d,P_C,E_C,w,N) - C_d
        EQ_C_f = E5_C_f(gamma_f,p_f,P_C,E_C,w,N) - C_f
        EQ_L_S = E6_L_S(N) - L
        EQ_Y   = E7_Y(C_d,X) - Y
        EQ_X   = E8_X(phi,p_d,p_f,E_X) - X
        EQ_p_f = E9_p_f()-p_f
        
        ## Returning a list of equations
        return [EQ_L_D,EQ_p_d,EQ_wN,EQ_C_d,EQ_C_f,EQ_L_S,EQ_Y,EQ_X,EQ_p_f]
   
    ## Solving the model  using the inital equilibrium as starting values fore the solver
    solution = optimize.root(CGEmodel,(ini_list[0],ini_list[1],ini_list[2],ini_list[3],ini_list[4],ini_list[5],ini_list[6],ini_list[7],ini_list[8]))
    
    ## Prints the status of the solver
    if status=='yes' :
        print(solution.message,"Success =",solution.success)
    
    ## Returning solution
    return solution

## Solving the model
zeroshock = CGEsolve(N,mu,gamma_d,gamma_f,phi,E_Y,E_C,E_X,ini_list,status='yes')
```

    The solution converged. Success = True



```python
## Function creating a dataframe with results
def CGEresults(solution,name,nice='yes') :
    """
    This function takes the results from CGEsolve function and stores the variable values in a data frame.
    
    solution = a result from the CGEsolve function
    name = name of the shock
    nice = yes if results are for printing in a table.
    """
    ## Saving the results in a dictionary - strings for presenting results
    if nice == 'yes' :
        results = [{"Shock" : name, 
                "$Y$"   : "%0.2f" % solution.x[0],
                "$L$"   : "%1.f"  % solution.x[1],
                "$w$"   : "%0.3f" % solution.x[2],
                "$C_d$" : "%0.2f" % solution.x[3],
                "$p_d$" : "%0.3f" % solution.x[5],
                "$C_f$" : "%0.2f" % solution.x[4],
                "$p_f$" : "%0.3f" % solution.x[6],
                "$P_C$" : "%0.3f" % solution.x[7],
                "$X$"   : "%0.2f" % solution.x[8]}]
    ## floats for figures 
    else :
        results = [{"Shock" : name, 
                "$Y$"   : solution.x[0],
                "$L$"   : solution.x[1],
                "$w$"   : solution.x[2],
                "$C_d$" : solution.x[3],
                "$p_d$" : solution.x[5],
                "$C_f$" : solution.x[4],
                "$p_f$" : solution.x[6],
                "$P_C$" : solution.x[7],
                "$X$"   : solution.x[8]}]
    ## Converting the results to a dataframe
    results = pd.DataFrame(results)
    
    ## Rearranging the columns
    cols = ['Shock','$L$','$Y$','$X$','$C_d$','$C_f$','$w$','$p_d$','$p_f$','$P_C$']
    results = results[cols]
    
    ## Returning the dataframe
    return results

## Printing the results
ini_eq.append(CGEresults(zeroshock,'Zero-shock',nice='yes'))
```




<div>
<style scoped>
    .dataframe tbody tr th:only-of-type {
        vertical-align: middle;
    }

    .dataframe tbody tr th {
        vertical-align: top;
    }

    .dataframe thead th {
        text-align: right;
    }
</style>
<table border="1" class="dataframe">
  <thead>
    <tr style="text-align: right;">
      <th></th>
      <th>Shock</th>
      <th>$L$</th>
      <th>$Y$</th>
      <th>$X$</th>
      <th>$C_d$</th>
      <th>$C_f$</th>
      <th>$w$</th>
      <th>$p_d$</th>
      <th>$p_f$</th>
      <th>$P_C$</th>
    </tr>
  </thead>
  <tbody>
    <tr>
      <th>0</th>
      <td>Initial eq.</td>
      <td>1000</td>
      <td>1000.00</td>
      <td>200.00</td>
      <td>800.00</td>
      <td>200.00</td>
      <td>1.000</td>
      <td>1.000</td>
      <td>1.000</td>
      <td>1.000</td>
    </tr>
    <tr>
      <th>0</th>
      <td>Zero-shock</td>
      <td>1000</td>
      <td>1000.00</td>
      <td>200.00</td>
      <td>800.00</td>
      <td>200.00</td>
      <td>1.000</td>
      <td>1.000</td>
      <td>1.000</td>
      <td>1.000</td>
    </tr>
  </tbody>
</table>
</div>



We have successfully solved the model. Our 'zero'-shock replicates the intial equilibrium. The model is correctly calibrated.

## The effects of a productivity shock and the price sensitivity of exports
In this section we examine how a 10 percent increase in productivity affects our model. The shok is implemented by increasing $\mu$ by 10 percent.


```python
prodshock = CGEsolve(N,1.1*mu,gamma_d,gamma_f,phi,E_Y,E_C,E_X,ini_list,status='yes')
```

    The solution converged. Success = True



```python
ini_eq.append(CGEresults(prodshock,'Prod-shock',nice='yes'))
```




<div>
<style scoped>
    .dataframe tbody tr th:only-of-type {
        vertical-align: middle;
    }

    .dataframe tbody tr th {
        vertical-align: top;
    }

    .dataframe thead th {
        text-align: right;
    }
</style>
<table border="1" class="dataframe">
  <thead>
    <tr style="text-align: right;">
      <th></th>
      <th>Shock</th>
      <th>$L$</th>
      <th>$Y$</th>
      <th>$X$</th>
      <th>$C_d$</th>
      <th>$C_f$</th>
      <th>$w$</th>
      <th>$p_d$</th>
      <th>$p_f$</th>
      <th>$P_C$</th>
    </tr>
  </thead>
  <tbody>
    <tr>
      <th>0</th>
      <td>Initial eq.</td>
      <td>1000</td>
      <td>1000.00</td>
      <td>200.00</td>
      <td>800.00</td>
      <td>200.00</td>
      <td>1.000</td>
      <td>1.000</td>
      <td>1.000</td>
      <td>1.000</td>
    </tr>
    <tr>
      <th>0</th>
      <td>Prod-shock</td>
      <td>1000</td>
      <td>1100.00</td>
      <td>221.83</td>
      <td>878.17</td>
      <td>217.28</td>
      <td>1.077</td>
      <td>0.979</td>
      <td>1.000</td>
      <td>0.984</td>
    </tr>
  </tbody>
</table>
</div>



The shock went as expected. A 10 percent increase in the productivity of the only input in a CES-function with constant returns to scale equals a 10 percent increase in domestic production. An increase in production also increases the export, as a larger supply lowers the price of the domestic good and thus increases exports. Also, a higher productivity leads to higher wages and thus higher consumption of both domestic and foregin goods. 

But how does the effect of the shock depend on $E_X$ in equation (8)? An inelastic (low values of $E_X$) foregin demand for domestic goods imply lower prices for domestic goods in order to get the market for domestic goods to clear. To see how dominant this effect is, we solve the model for a range of $E_X$'s. 


```python
## Solving the models for E_X=2
shock = CGEsolve(N,1.1*mu,gamma_d,gamma_f,phi,E_Y,E_C,2,ini_list,status='yes')
results = CGEresults(shock,2,nice = 'no')
```

    The solution converged. Success = True



```python
## Solving the model for E_X=2-30 - [we don't print the status of the solver, but we have checked. 
##                                   The solver behaves nicely for all values of E_X]
for i in range(21,301) :
    #print(0.1*i)
    shock = CGEsolve(N,1.1*mu,gamma_d,gamma_f,phi,E_Y,E_C,0.1*i,ini_list,status='no')
    results = results.append(CGEresults(shock,0.1*i,nice='no'))
```

We present the results in an interactive figure.


```python
varbls = ['$L$','$Y$','$X$','$C_d$','$C_f$','$w$','$p_d$','$p_f$','$P_C$']

def fig_plt(var) :
    y = results[var].values
    x = results['Shock'].values
    plt.plot(x,y)
    plt.yticks(np.arange(0.95*y.min(), 1.05*y.max(), (1.05*y.max()-0.95*y.min())/4))
    return plt.show()

var_select = wd.Dropdown(options=varbls, description='Variables')
wd.interact(fig_plt,var=var_select);
```


    interactive(children=(Dropdown(description='Variables', options=('$L$', '$Y$', '$X$', '$C_d$', '$C_f$', '$w$',…


The most interesting figures are those of $X$ and $p_d$. It shows that an inelastic foregin demand implies, that the domestic producers have to lower their prices by a lot to clear the goods market. The lower prices affects the wages and thus the consumption of both foregin and domestic goods in equilibria with different $E_X$.

237657f6af57442af5cbef59d9339d1fc3fed83f
ipynb
Jupyter Notebook
modelproject/CGE_model.ipynb
NumEconCopenhagen/projects-2019-faetter-br
80c5afb7496dcd756f12cde9c0d988a891cf9619
MIT
modelproject/CGE_model.ipynb
NumEconCopenhagen/projects-2019-faetter-br
80c5afb7496dcd756f12cde9c0d988a891cf9619
MIT
2019-04-10T10:38:49.000Z
2019-05-14T19:22:39.000Z
modelproject/CGE_model.ipynb
NumEconCopenhagen/projects-2019-faetter-br
80c5afb7496dcd756f12cde9c0d988a891cf9619
MIT
2019-04-02T08:27:23.000Z
2019-04-02T08:27:23.000Z
Qwen/Qwen-72B
1. YES
2. YES