Results

In this page we report the results from running the package with the official scenarios and compare the output of DICEModel.jl with the official GAMS output as given in the DICE 2023 Introduction and User's Manual (v3.1.2, May 15, 2024)

The scenarios considered are:

cbopt: The C/B optimal scenario
t2c: The temperature constrained to max 2 °C scenario
t15c: The temperature constrained to max 1.5 °C scenario
altdam: The alternative damage scenario
parisext: The Paris extended scenario
base: The base (current policies) scenario
r5: The scenario with discount rate of 5%
r4: The scenario with discount rate of 4%
r3: The scenario with discount rate of 3%
r2: The scenario with discount rate of 2%
r1: The scenario with discount rate of 1%

Results

We start by loading the required packages.

Show code

using Pkg
Pkg.activate(joinpath(@__DIR__,".."))
#Pkg.add(["DICEModel", "CSV","DataFrames","Plots"]) # Run once and comment this back
using DICEModel, CSV, DataFrames, Plots, XLSX

  Activating project at `~/work/DICEModel.jl/DICEModel.jl/docs`

We can now define the scenarios to run and their specific graphical attributed when plotted. Finally we call the function run_dice_scenario with each of them and save the output in the results dictionary (keyed by the scenario name):

Show code

# Scenarios to run:
scenarios = ["cbopt","t2c","t15c","altdam","parisext","base","r5","r4","r3","r2","r1"]
plot_attributes = Dict(
    "cbopt"    => (label="C/B optimal", linestyle=:solid, colour=:red, markershape=:circle, markercolor=:white, markerstrokecolor=:red),
    "t2c"      => (label="T < 2 °C", linestyle=:dash, colour=:green, markershape=:cross, markercolor=:green, markerstrokecolor=:green),
    "t15c"     => (label="T < 1.5 °C", linestyle=:dot, colour=:lime, markershape=:dtriangle, markercolor=:lime, markerstrokecolor=:lime),
    "altdam"   => (label="Alt. damage", linestyle=:solid, colour=:darkgoldenrod4, markershape=:none, markercolor=:darkgoldenrod4, markerstrokecolor=:darkgoldenrod4),
    "parisext" => (label="Paris ext", linestyle=:dash, colour=:darkorange, markershape=:cross, markercolor=:darkorange, markerstrokecolor=:darkorange),
    "base"     => (label="Base", linestyle=:solid, colour=:yellow, markershape=:utriangle, markercolor=:yellow, markerstrokecolor=:yellow),
    "r5"       => (label="R = 5%", linestyle=:solid, colour=:darkblue, markershape=:diamond, markercolor=:darkblue, markerstrokecolor=:darkblue),
    "r4"       => (label="R = 4%", linestyle=:dash, colour=:blue, markershape=:diamond, markercolor=:blue, markerstrokecolor=:blue),
    "r3"       => (label="R = 3%", linestyle=:solid, colour=:deepskyblue4, markershape=:diamond, markercolor=:deepskyblue4, markerstrokecolor=:deepskyblue4),
    "r2"       => (label="R = 2%", linestyle=:dash, colour=:steelblue, markershape=:diamond, markercolor=:steelblue, markerstrokecolor=:steelblue),
    "r1"       => (label="R = 1%", linestyle=:solid, colour=:lightskyblue, markershape=:diamond, markercolor=:lightskyblue, markerstrokecolor=:lightskyblue)
)
# Run the scenarios and collect their results in the "results" dictionary
results = Dict([s => run_dice_scenario(s) for s in scenarios])
times  = results["cbopt"].times.+2020


******************************************************************************
This program contains Ipopt, a library for large-scale nonlinear optimization.
 Ipopt is released as open source code under the Eclipse Public License (EPL).
         For more information visit https://github.com/coin-or/Ipopt
******************************************************************************

Main results and comparison with the official GAMS version

For the variables ECO2, ppm, TATM, MIU, CPRICE and scc, we show the official values, the ones computed with DICEModel.jl and we plot these values.

`ECO2` : CO₂ emissions [GtCO₂/yr]

Official GAMS outputs:

Show code

# Data from the DICE2023 manual
ECO2_gams = CSV.read(IOBuffer("""
scen      2020  2025  2050  2100
cbopt     42.9  42.9  37.1  15.9
t2c       42.9  42.9  27.2   1.2
t15c      42.9  13.1   5.7   0.0
altdam    42.9  42.7  20.9   0.0
parisext  42.9  43.3  44.4  42.3
base      42.9  44.9  54.6  75.7
r5        42.8  42.5  42.2  37.6
r4        44.1  43.9  39.3  28.9
r3        45.6  45.3  33.5  15.4
r2        46.8  46.7  22.5   0.0
r1        46.8  46.9  19.2   0.0
"""), DataFrame, delim=" ", ignorerepeated=true)

11×5 DataFrame

Row	scen	2020	2025	2050	2100
	String15	Float64	Float64	Float64	Float64
1	cbopt	42.9	42.9	37.1	15.9
2	t2c	42.9	42.9	27.2	1.2
3	t15c	42.9	13.1	5.7	0.0
4	altdam	42.9	42.7	20.9	0.0
5	parisext	42.9	43.3	44.4	42.3
6	base	42.9	44.9	54.6	75.7
7	r5	42.8	42.5	42.2	37.6
8	r4	44.1	43.9	39.3	28.9
9	r3	45.6	45.3	33.5	15.4
10	r2	46.8	46.7	22.5	0.0
11	r1	46.8	46.9	19.2	0.0

DICEModel.jl output

Show code

ECO2 = DataFrame([
    scenarios,
    [results[s].ECO2[1] for s in scenarios],
    [results[s].ECO2[2] for s in scenarios],
    [results[s].ECO2[7] for s in scenarios],
    [results[s].ECO2[17] for s in scenarios],
    ],
    ["scen","2020","2025","2050","2100"]
)

11×5 DataFrame

Row	scen	2020	2025	2050	2100
	String	Float64	Float64	Float64	Float64
1	cbopt	42.9412	42.9103	37.0559	15.892
2	t2c	42.9412	42.8565	27.1827	1.20346
3	t15c	42.9412	13.06	5.66113	-7.15383e-7
4	altdam	42.9412	42.6845	20.8841	-9.29042e-7
5	parisext	42.9412	43.349	44.3675	42.3349
6	base	42.8829	44.8858	54.6281	75.6928
7	r5	42.7502	42.4828	42.2009	37.7465
8	r4	44.0774	43.8959	39.3474	28.885
9	r3	45.5667	45.348	33.5202	15.4198
10	r2	46.7863	46.6775	22.4688	-1.0527e-6
11	r1	47.1155	46.8666	19.1864	-1.11658e-6

Differences:

Show code

ECO2_diff = copy(ECO2_gams);
ECO2_diff[:,2:end] .= ECO2_gams[:,2:end] .- ECO2[:,2:end]

11×5 DataFrame

Row	scen	2020	2025	2050	2100
	String15	Float64	Float64	Float64	Float64
1	cbopt	-0.0411997	-0.010322	0.0440893	0.00799729
2	t2c	-0.0411997	0.0435293	0.0173196	-0.00346045
3	t15c	-0.0411996	0.0399781	0.0388712	7.15383e-7
4	altdam	-0.0411997	0.0155168	0.0159117	9.29042e-7
5	parisext	-0.0411997	-0.0490415	0.0325285	-0.0348967
6	base	0.0171207	0.0142431	-0.0280874	0.00719638
7	r5	0.0497826	0.0172372	-0.00086204	-0.146482
8	r4	0.0226438	0.00406249	-0.0473563	0.0150089
9	r3	0.0332938	-0.0480203	-0.0202132	-0.019788
10	r2	0.0136577	0.0224956	0.0311961	1.0527e-6
11	r1	-0.315531	0.0334078	0.0136131	1.11658e-6

There are only two minor differences in r5 for 2100 and, above all, for r1 in 2020.

Plot:

Show code

scenarios_plot=["base","cbopt","parisext","t2c","r2","r1"]
p = plot()
for(i,s) in enumerate(scenarios_plot)
    if i == 1
        global p = plot(times[1:11] ,results[s].ECO2[1:11],ylim=(0,80), title="CO₂ emissions",ylabel="GtCO₂/yr";plot_attributes[s]...);
    else
        plot!(times[1:11] ,results[s].ECO2[1:11];plot_attributes[s]...)
    end
end

`ppm`: CO₂ concentration [ppm]

Official GAMS outputs:

Show code

# Data from the DICE2023 manual
ppm_gams = CSV.read(IOBuffer("""
scen       2020  2025  2050  2100   2150
cbopt     416.2 429.9 487.8 569.2  497.9
t2c       416.2 429.9 474.7 474.7  437.9
altdam    416.2 429.8 466.7 458.5  401.0
parisext  416.2 430.1 501.3 652.5  763.5
base      416.2 430.9 517.7 774.9 1144.0
r5        416.2 429.7 495.7 635.5  671.6
r4        416.2 430.4 491.8 605.3  592.0
r3        416.2 431.2 484.0 555.1  494.2
r2        416.2 431.9 473.8 484.9  419.3
r1        416.2 432.0 473.3 449.5  389.3
"""), DataFrame, delim=" ", ignorerepeated=true)

10×6 DataFrame

Row	scen	2020	2025	2050	2100	2150
	String15	Float64	Float64	Float64	Float64	Float64
1	cbopt	416.2	429.9	487.8	569.2	497.9
2	t2c	416.2	429.9	474.7	474.7	437.9
3	altdam	416.2	429.8	466.7	458.5	401.0
4	parisext	416.2	430.1	501.3	652.5	763.5
5	base	416.2	430.9	517.7	774.9	1144.0
6	r5	416.2	429.7	495.7	635.5	671.6
7	r4	416.2	430.4	491.8	605.3	592.0
8	r3	416.2	431.2	484.0	555.1	494.2
9	r2	416.2	431.9	473.8	484.9	419.3
10	r1	416.2	432.0	473.3	449.5	389.3

DICEModel.jl output

Show code

scenarios_ppm = ["cbopt","t2c","altdam","parisext","base","r5","r4","r3","r2","r1"]
ppm = DataFrame([
    scenarios_ppm,
    [results[s].ppm[1] for s in scenarios_ppm],
    [results[s].ppm[2] for s in scenarios_ppm],
    [results[s].ppm[7] for s in scenarios_ppm],
    [results[s].ppm[17] for s in scenarios_ppm],
    [results[s].ppm[27] for s in scenarios_ppm],
    ],
    ["scen","2020","2025","2050","2100","2150"]
)

10×6 DataFrame

Row	scen	2020	2025	2050	2100	2150
	String	Float64	Float64	Float64	Float64	Float64
1	cbopt	416.203	429.915	487.795	569.243	497.932
2	t2c	416.203	429.888	474.688	474.749	437.865
3	altdam	416.203	429.802	466.698	458.507	401.008
4	parisext	416.203	430.134	501.294	652.492	763.456
5	base	416.203	430.898	517.714	774.93	1143.99
6	r5	416.203	429.697	495.734	635.758	673.822
7	r4	416.203	430.433	491.782	605.346	591.963
8	r3	416.203	431.192	484.031	555.081	494.245
9	r2	416.203	431.884	473.81	484.93	419.344
10	r1	416.203	431.986	473.289	449.544	389.345

Differences:

Show code

ppm_diff = copy(ppm_gams);
ppm_diff[:,2:end] .= ppm_gams[:,2:end] .- ppm[:,2:end]

10×6 DataFrame

Row	scen	2020	2025	2050	2100	2150
	String15	Float64	Float64	Float64	Float64	Float64
1	cbopt	-0.00319315	-0.014811	0.00507806	-0.0425309	-0.0324887
2	t2c	-0.00319315	0.0120542	0.0116026	-0.0488929	0.0351935
3	altdam	-0.00319315	-0.00214141	0.00183298	-0.00691756	-0.00801501
4	parisext	-0.00319315	-0.0336583	0.00621909	0.00768185	0.0439192
5	base	-0.00319315	0.00156887	-0.0142669	-0.0296323	0.00881173
6	r5	-0.00319315	0.00286689	-0.0343345	-0.257502	-2.22231
7	r4	-0.00319315	-0.032674	0.0177718	-0.0458003	0.0366566
8	r3	-0.00319315	0.00825996	-0.0306922	0.0188462	-0.0451211
9	r2	-0.00319315	0.0163165	-0.00950111	-0.0300991	-0.0435679
10	r1	-0.00319315	0.0142086	0.010667	-0.0440142	-0.0447197

Again, the only minor difference is for r5 in 2150.

Plot:

Show code

scenarios_plot=["base","cbopt","parisext","t2c","r2","r1"]
for(i,s) in enumerate(scenarios_plot)
    if i == 1
        global p = plot(times[1:17] ,results[s].ppm[1:17],ylim=(400,800), title="CO₂ concentrations",ylabel="ppm";plot_attributes[s]...);
    else
        plot!(times[1:17] ,results[s].ppm[1:17];plot_attributes[s]...)
    end
end

`TATM`: Global Temperatures [°C diff since 1765]

Official GAMS outputs:

Show code

# Data from the DICE2023 manual
TATM_gams = CSV.read(IOBuffer("""
scen       2020  2025  2050  2100  2150
cbopt      1.25  1.42  1.92  2.58  2.29
t2c        1.25  1.42  1.85  2.00  1.86
altdam     1.25  1.42  1.81  1.89  1.58
parisext   1.25  1.43  2.01  3.00  3.61
base       1.25  1.43  2.10  3.55  4.91
r5         1.25  1.42  1.97  2.93  3.24
r4         1.25  1.43  1.95  2.77  2.84
r3         1.25  1.43  1.90  2.49  2.26
r2         1.25  1.43  1.84  2.07  1.73
r1         1.25  1.43  1.84  1.81  1.49
"""), DataFrame, delim=" ", ignorerepeated=true)

10×6 DataFrame

Row	scen	2020	2025	2050	2100	2150
	String15	Float64	Float64	Float64	Float64	Float64
1	cbopt	1.25	1.42	1.92	2.58	2.29
2	t2c	1.25	1.42	1.85	2.0	1.86
3	altdam	1.25	1.42	1.81	1.89	1.58
4	parisext	1.25	1.43	2.01	3.0	3.61
5	base	1.25	1.43	2.1	3.55	4.91
6	r5	1.25	1.42	1.97	2.93	3.24
7	r4	1.25	1.43	1.95	2.77	2.84
8	r3	1.25	1.43	1.9	2.49	2.26
9	r2	1.25	1.43	1.84	2.07	1.73
10	r1	1.25	1.43	1.84	1.81	1.49

DICEModel.jl output

Show code

scenarios_TATM = ["cbopt","t2c","altdam","parisext","base","r5","r4","r3","r2","r1"]
TATM = DataFrame([
    scenarios_TATM,
    [results[s].TATM[1] for s in scenarios_TATM],
    [results[s].TATM[2] for s in scenarios_TATM],
    [results[s].TATM[7] for s in scenarios_TATM],
    [results[s].TATM[17] for s in scenarios_TATM],
    [results[s].TATM[27] for s in scenarios_TATM],
    ],
    ["scen","2020","2025","2050","2100","2150"]
)

10×6 DataFrame

Row	scen	2020	2025	2050	2100	2150
	String	Float64	Float64	Float64	Float64	Float64
1	cbopt	1.24715	1.42484	1.92441	2.58054	2.28903
2	t2c	1.24715	1.42473	1.85038	2.0	1.85951
3	altdam	1.24715	1.42437	1.80619	1.88723	1.57725
4	parisext	1.24715	1.42576	2.00544	3.00231	3.60536
5	base	1.24715	1.42893	2.09845	3.55492	4.91077
6	r5	1.24715	1.42393	1.97261	2.92756	3.25007
7	r4	1.24715	1.42701	1.94585	2.76606	2.84041
8	r3	1.24715	1.43019	1.89857	2.48968	2.25982
9	r2	1.24715	1.43308	1.84487	2.06604	1.72796
10	r1	1.24715	1.43351	1.84504	1.81495	1.48749

Differences:

Show code

TATM_diff = copy(TATM_gams);
TATM_diff[:,2:end] .= TATM_gams[:,2:end] .- TATM[:,2:end]

10×6 DataFrame

Row	scen	2020	2025	2050	2100	2150
	String15	Float64	Float64	Float64	Float64	Float64
1	cbopt	0.002846	-0.00484409	-0.00440847	-0.00054271	0.000974924
2	t2c	0.002846	-0.00473142	-0.000380345	-1.95338e-8	0.000493873
3	altdam	0.002846	-0.00437152	0.00380663	0.00277346	0.00274842
4	parisext	0.002846	0.00423836	0.00456118	-0.00230776	0.00464193
5	base	0.002846	0.00107463	0.00154796	-0.00491501	-0.000773599
6	r5	0.002846	-0.00393098	-0.00261055	0.0024408	-0.0100711
7	r4	0.002846	0.00298545	0.00415107	0.00394109	-0.000413549
8	r3	0.002846	-0.000191228	0.00142958	0.000324026	0.000182754
9	r2	0.002846	-0.00308213	-0.00487368	0.0039635	0.00204289
10	r1	0.002846	-0.00350834	-0.00504228	-0.00494626	0.00251218

Again, the only minor difference is for r5 in 2150 (3.25 instead of 3.24)

Plot:

Show code

scenarios_plot=["base","cbopt","parisext","t2c","r4","r2"]
for(i,s) in enumerate(scenarios_plot)
    if i == 1
        global p = plot(times[1:17] ,results[s].TATM[1:17],ylim=(0.0,4.0), title="Global Temperatures",ylabel="°C from 1765";plot_attributes[s]...);
    else
        plot!(times[1:17] ,results[s].TATM[1:17];plot_attributes[s]...)
    end
end

`MIU`: Emission control rate [%]

Official GAMS outputs:

Show code

# Data from the DICE2023 manual
MIU_gams = CSV.read(IOBuffer("""
scen    2020 2030 2040 2050 2060 2100
cbopt      5   24   31   39   46   84
t2c        5   24   42   55   69   99
altdam     5   24   48   65   76  100
parisext   5   13   21   27   33   57
base       5    6    8   10   12   22
r5         5   19   23   29   34   60
r4         5   24   29   36   42   70
r3         5   24   39   47   54   85
r2         5   24   48   66   73  100
r1         5   24   48   72   90  100
"""), DataFrame, delim=" ", ignorerepeated=true)

10×7 DataFrame

Row	scen	2020	2030	2040	2050	2060	2100
	String15	Int64	Int64	Int64	Int64	Int64	Int64
1	cbopt	5	24	31	39	46	84
2	t2c	5	24	42	55	69	99
3	altdam	5	24	48	65	76	100
4	parisext	5	13	21	27	33	57
5	base	5	6	8	10	12	22
6	r5	5	19	23	29	34	60
7	r4	5	24	29	36	42	70
8	r3	5	24	39	47	54	85
9	r2	5	24	48	66	73	100
10	r1	5	24	48	72	90	100

DICEModel.jl output

Show code

scenarios_MIU = ["cbopt","t2c","altdam","parisext","base","r5","r4","r3","r2","r1"]
MIU = DataFrame([
    scenarios_MIU,
    [results[s].MIU[1] for s in scenarios_MIU] .* 100,
    [results[s].MIU[3] for s in scenarios_MIU] .* 100,
    [results[s].MIU[5] for s in scenarios_MIU] .* 100,
    [results[s].MIU[7] for s in scenarios_MIU] .* 100,
    [results[s].MIU[9] for s in scenarios_MIU] .* 100,
    [results[s].MIU[17] for s in scenarios_MIU] .* 100,
    ],
    ["scen","2020","2030","2040","2050","2060","2100"]
)

10×7 DataFrame

Row	scen	2020	2030	2040	2050	2060	2100
	String	Float64	Float64	Float64	Float64	Float64	Float64
1	cbopt	5.0	24.0	30.9921	39.1944	46.45	84.0037
2	t2c	5.0	24.0	41.801	55.2869	69.2618	98.7874
3	altdam	5.0	24.0	48.0	65.3108	75.6007	100.0
4	parisext	5.0	13.0	21.0	27.0	33.0	57.0
5	base	5.12902	6.37092	7.91352	9.82964	11.5419	21.9397
6	r5	5.0	18.4969	23.2351	28.9657	33.895	59.8339
7	r4	5.0	23.6722	29.2991	36.0564	41.7016	70.4212
8	r3	5.0	24.0	39.3833	47.4538	53.8342	84.848
9	r2	5.0	24.0	48.0	66.045	72.5916	100.0
10	r1	5.0	24.0	48.0	72.0	90.0	100.0

Differences:

Show code

MIU_diff = copy(MIU);
MIU_diff[:,2:end] .= MIU_gams[:,2:end] .- MIU[:,2:end]

10×7 DataFrame

Row	scen	2020	2030	2040	2050	2060	2100
	String	Float64	Float64	Float64	Float64	Float64	Float64
1	cbopt	-8.70517e-7	2.86129e-7	0.00793158	-0.194426	-0.449981	-0.00367471
2	t2c	-9.20874e-7	-9.55094e-7	0.199016	-0.286948	-0.261819	0.212592
3	altdam	-9.58301e-7	-9.83517e-7	-9.27644e-7	-0.310777	0.399331	-9.45696e-7
4	parisext	-9.23627e-7	-9.70824e-7	-9.63689e-7	-9.62881e-7	-9.61297e-7	-9.63202e-7
5	base	-0.129025	-0.370923	0.0864759	0.170363	0.458135	0.0603352
6	r5	-9.77287e-7	0.503117	-0.235148	0.0342768	0.104978	0.16612
7	r4	-9.90236e-7	0.327775	-0.299091	-0.0563595	0.298353	-0.421161
8	r3	-9.95254e-7	-9.97168e-7	-0.383298	-0.453829	0.165817	0.152027
9	r2	-9.97388e-7	-9.99252e-7	-9.97949e-7	-0.0450235	0.408357	-9.88648e-7
10	r1	-9.98526e-7	-9.99798e-7	-9.99762e-7	-9.99692e-7	-9.99568e-7	-9.99738e-7

No differences here !

Plot:

Show code

scenarios_plot=["cbopt","parisext","t2c","altdam"]
for(i,s) in enumerate(scenarios_plot)
    if i == 1
        global p = plot(times[1:17] ,results[s].MIU[1:17] .* 100,ylim=(0.0,100), title="Emission control rate",ylabel="%";plot_attributes[s]...);
    else
        plot!(times[1:17] ,results[s].MIU[1:17] .* 100;plot_attributes[s]...)
    end
end

`CPRICE`: Carbon price [2019$ / tCO₂]

Plot:

Show code

scenarios_plot=["cbopt","parisext","t2c","r3"]
for(i,s) in enumerate(scenarios_plot)
    if i == 1
        global p = plot(times[1:9] ,results[s].CPRICE[1:9],ylim=(0,350), title="Carbon price",ylabel="2019\$ / tCO₂";plot_attributes[s]...);
    else
        plot!(times[1:9] ,results[s].CPRICE[1:9];plot_attributes[s]...)
    end
end

Official GAMS outputs:

Show code

# Data from the DICE2023 manual
scc_gams = CSV.read(IOBuffer("""
scen       2020  2025  2050
cbopt        50    59   125
t2c          75    89   213
t15c       3557  4185 16552
altdam      124   146   281
parisext     61    72   159
base         66    78   175
r5           32    37    74
r4           49    58   107
r3           87   102   172
r2          176   207   302
r1          485   571   695
"""), DataFrame, delim=" ", ignorerepeated=true)

11×4 DataFrame

Row	scen	2020	2025	2050
	String15	Int64	Int64	Int64
1	cbopt	50	59	125
2	t2c	75	89	213
3	t15c	3557	4185	16552
4	altdam	124	146	281
5	parisext	61	72	159
6	base	66	78	175
7	r5	32	37	74
8	r4	49	58	107
9	r3	87	102	172
10	r2	176	207	302
11	r1	485	571	695

DICEModel.jl output

Show code

scenarios_scc = ["cbopt","t2c","t15c", "altdam","parisext","base","r5","r4","r3","r2","r1"]
scc = DataFrame([
    scenarios_scc,
    [results[s].scc[1] for s in scenarios_scc],
    [results[s].scc[2] for s in scenarios_scc],
    [results[s].scc[7] for s in scenarios_scc],
    ],
    ["scen","2020","2025","2050"]
)

11×4 DataFrame

Row	scen	2020	2025	2050
	String	Float64	Float64	Float64
1	cbopt	9.21936	59.0669	124.542
2	t2c	13.1877	89.034	213.132
3	t15c	1773.01	4184.56	16551.8
4	altdam	17.8033	145.918	281.084
5	parisext	15.9524	71.5199	159.426
6	base	21.3915	77.6798	175.028
7	r5	5.64394	37.2822	73.5004
8	r4	10.1477	57.685	106.929
9	r3	18.5177	101.816	171.605
10	r2	30.5636	206.997	302.131
11	r1	46.4243	570.813	695.261

Differences:

Show code

scc_diff = copy(scc);
scc_diff[:,2:end] .= scc_gams[:,2:end] .- scc[:,2:end]

11×4 DataFrame

Row	scen	2020	2025	2050
	String	Float64	Float64	Float64
1	cbopt	40.7806	-0.0669242	0.457849
2	t2c	61.8123	-0.034009	-0.131736
3	t15c	1783.99	0.444098	0.185102
4	altdam	106.197	0.0819754	-0.0837618
5	parisext	45.0476	0.480074	-0.425736
6	base	44.6085	0.320193	-0.0284349
7	r5	26.3561	-0.282159	0.499588
8	r4	38.8523	0.315	0.0711168
9	r3	68.4823	0.18418	0.394926
10	r2	145.436	0.0027403	-0.130938
11	r1	438.576	0.187235	-0.26097

This is the only part that still needs to be checked, as there are significant differences for the base year (2020). For the other years DICEModel.jl provides identical results (up to the approximation of the data published) than the official GAMS version.

Plot:

Show code

scenarios_plot=["base","cbopt","t2c","altdam","r4","r2"]
for(i,s) in enumerate(scenarios_plot)
    if i == 1
        global p = plot(times[1:7] ,results[s].scc[1:7],ylim=(0,400), title="Social cost of carbon",ylabel="2019\$ / tCO₂";plot_attributes[s]...);
    else
        plot!(times[1:7] ,results[s].scc[1:7];plot_attributes[s]...)
    end
end

Detailed model output

The whole model output (variables and post-processing computed values) can be retrieved in the Excel file below:

Show code

keys(results["cbopt"])
out_vars = [v for v in keys(results["cbopt"]) if typeof(results["cbopt"][v]) <: Vector{<: Number}]

all_results = DataFrame(
    scenario = String[],
    variable = String[],
    year     = Int64[],
    value    = Float64[],
)

for s in scenarios, v in out_vars, ti in 1:length(times)
    push!(all_results,[s,string(v),times[ti],results[s][v][ti]])
end
all_results = unstack(all_results,"year","value")

XLSX.writetable("DICEModelDetailedResults.xlsx",
all_results=(collect(DataFrames.eachcol(all_results)),DataFrames.names(all_results))
)

DICEModelDetailedResults.xlsx