Model close to Miranda and Glauber (1995). Countries $a$ and $b$ are indicated by the subscript $i \in \left\{a,b\right\}$. When variables corresponding to both countries appear in the same equation, the foreign country is indicated by the subscript $j$.

## Model's structure

Response variables Storage ($S_{i}$), Planned production ($H_{i}$), Price ($P_{i}$) and Export ($X_{i}$).

State variable Availability ($A_{i}$).

Shock Productivity shocks ($\epsilon_{i}$).

Parameters Unit physical storage cost ($k$), Depreciation share ($\delta$), Interest rate ($r$), Scale parameter for the production cost function ($h$), Inverse of supply elasticity ($\mu$), Demand price elasticity ($\alpha$), Scale parameter for demand function ($\gamma_i$), Trade cost ($\theta$) and Tariff ($\tau_{i}$).

Equilibrium equations

For $i \in \left\{a,b\right\}$ and $j \ne i$:

$S_{it}: S_{it}\ge 0 \quad \perp \quad P_{it}+k-\frac{1-\delta}{1+r}\mathrm{E}_{t}\left(P_{it+1}\right)\ge 0,$

$H_{it}: \frac{1}{1+r}\mathrm{E}_{t}\left(P_{it+1}\epsilon_{it+1}\right)=h {H_{it}}^{\mu},$

$P_{it}: A_{it}+X_{jt}=\gamma_i{P_{it}}^{\alpha}+S_{it}+X_{it},$

$X_{it}: X_{it}\ge 0 \quad \perp \quad P_{it}+\theta+\tau_{j}\ge P_{jt}.$

Transition equation

For $i \in \left\{a,b\right\}$:

$A_{it}: A_{it}=\left(1-\delta\right)S_{it-1}+H_{it-1}\epsilon_{it}.$

## Writing the model

The model is defined in a Yaml file: sto6.yaml.

## Create the model object

Mu                = [1 1];
sigma             = [0.05 0;
0    0.05];
model = recsmodel('sto6.yaml',struct('Mu',Mu,'Sigma',sigma^2,'order',7));

Deterministic steady state (different from first guess, max(|delta|)=0.31777)
State variables:
1.04	1.057

Response variables:
0	4.378e-20	1.04	1.057	1.218	1.318	0.07883	0

Expectations variables:
1.218	1.318	1.218	1.318



## Define approximation space

[interp,s] = recsinterpinit(15,0.73*model.sss,2*model.sss);


## Find a first guess through the perfect foresight solution

interp = recsFirstGuess(interp,model,s,model.sss,model.xss,struct('T',5));


## Solve for rational expectations

[interp,x] = recsSolveREE(interp,model,s);

Successive approximation
Major	 Minor	Residual
0	     0	7.73E-01 (Input point)
1	     1	2.14E-01
2	     1	9.75E-02
3	     1	7.31E-02
4	     1	4.34E-02
5	     1	1.69E-02
6	     1	4.49E-03
7	     1	1.00E-03
8	     1	2.10E-04
9	     1	4.32E-05
10	     1	8.82E-06
11	     1	1.80E-06
12	     1	3.74E-07
13	     1	8.46E-08
14	     1	0.00E+00
Solution found - Residual lower than absolute tolerance


## Simulate the model

[~,~,~,stat] = recsSimul(model,interp,model.sss(ones(1E4,1),:),100);

Statistics from simulated variables (excluding the first 20 observations):
Moments
Mean      Std. Dev. Skewness  Kurtosis  Min       Max       %LB       %UB
1.0615    0.0560    0.0792    3.0577    0.8200    1.3420
1.0588    0.0532    0.0022    3.0010    0.7835    1.3038

Columns 1 through 7

0.0214    0.0269    1.7753    6.5462         0    0.2633   20.2164
0.0014    0.0050    6.2398   57.8325         0    0.1138   70.0932
1.0406    0.0057   -1.4955    5.0697    1.0027    1.0463         0
1.0574    0.0055   -1.5053    5.0941    1.0214    1.0628         0
1.2221    0.1754    1.5256    6.0478    0.9065    3.1226         0
1.3241    0.1719    1.5480    6.5039    0.9988    3.2159         0
0.0779    0.0388    0.3320    2.6825         0    0.2682    0.1539
0.0002    0.0009   23.5274  799.8128         0    0.0620   44.8194

Column 8

0
0
0
0
0
0
0
0

Correlation
Columns 1 through 7

1.0000   -0.0327    0.6834   -0.0845   -0.6456   -0.6457   -0.6615
-0.0327    1.0000    0.5312    0.5269   -0.6049   -0.6035   -0.6360
0.6834    0.5312    1.0000    0.1960   -0.9817   -0.9822   -0.6589
-0.0845    0.5269    0.1960    1.0000   -0.3655   -0.3637   -0.2183
-0.6456   -0.6049   -0.9817   -0.3655    1.0000    1.0000    0.6842
-0.6457   -0.6035   -0.9822   -0.3637    1.0000    1.0000    0.6828
-0.6615   -0.6360   -0.6589   -0.2183    0.6842    0.6828    1.0000
-0.6695   -0.6355   -0.6734   -0.2212    0.6971    0.6958    0.9955
0.5226   -0.8496   -0.1820   -0.4211    0.2431    0.2423    0.1562
-0.2083    0.1136   -0.0607    0.0776    0.0556    0.0558    0.1083

Columns 8 through 10

-0.6695    0.5226   -0.2083
-0.6355   -0.8496    0.1136
-0.6734   -0.1820   -0.0607
-0.2212   -0.4211    0.0776
0.6971    0.2431    0.0556
0.6958    0.2423    0.0558
0.9955    0.1562    0.1083
1.0000    0.1500    0.0475
0.1500    1.0000   -0.1605
0.0475   -0.1605    1.0000

Autocorrelation
1         2         3         4         5
0.2296    0.0295   -0.0108   -0.0176   -0.0205
-0.0287   -0.0247   -0.0152   -0.0123   -0.0127

0.1978    0.0227   -0.0113   -0.0175   -0.0189
0.0059   -0.0161   -0.0147   -0.0138   -0.0121
0.1821    0.0198   -0.0113   -0.0176   -0.0182
0.1824    0.0198   -0.0113   -0.0176   -0.0182
0.1091    0.0078   -0.0132   -0.0166   -0.0165
0.1153    0.0091   -0.0131   -0.0167   -0.0165
-0.0304   -0.0263   -0.0150   -0.0115   -0.0133
-0.0055   -0.0118   -0.0132   -0.0128   -0.0141



## References

Miranda, M. J. and Glauber, J. W. (1995). Solving stochastic models of competitive storage and trade by Chebychev collocation methods. Agricultural and Resource Economics Review, 24(1), 70-77.