diff --git a/results/01_initial_experiment.txt b/results/01_initial_experiment.txt new file mode 100644 index 0000000..89996e8 --- /dev/null +++ b/results/01_initial_experiment.txt @@ -0,0 +1,43 @@ +parameters values +beta 0.528253602365166 +q10_rh 1.76 +diff 2.3843958670058 + +Objective Function Value: 2.3 +Counts: 15 +Counts: 15 +Convergence: 0 +Messages: CONVERGENCE: REL_REDUCTION_OF_F <= FACTR*EPSMCH + +CO2 MSE: 4.57 +T MSE: 0.0357 +RMSE: 2.14 + +***Key Metrics*** +TCRE: 1.51 +TCR: 1.77 + +***Historical Warming and ERF*** +GSAT Warming: 0.731 +Ocean Heat Content Change: 472 +Total Aerosol ERF: -1.24 +WMGHG ERF: 3.87 +Methane ERF: 0.54 + +***Future Warming*** +scenario start end GSAT +ssp119 2021 2040 0.73 +ssp119 2041 2060 0.895 +ssp119 2081 2100 0.717 +ssp126 2021 2040 0.745 +ssp126 2041 2060 1.08 +ssp126 2081 2100 1.1 +ssp245 2021 2040 0.75 +ssp245 2041 2060 1.29 +ssp245 2081 2100 1.98 +ssp370 2021 2040 0.761 +ssp370 2041 2060 1.43 +ssp370 2081 2100 2.94 +ssp585 2021 2040 0.877 +ssp585 2041 2060 1.74 +ssp585 2081 2100 3.79 diff --git a/results/02_smooth_temps_3yr.txt b/results/02_smooth_temps_3yr.txt new file mode 100644 index 0000000..a62d312 --- /dev/null +++ b/results/02_smooth_temps_3yr.txt @@ -0,0 +1,43 @@ +parameters values +beta 0.528263071532207 +q10_rh 1.76 +diff 2.37997014792955 + +Objective Function Value: 2.3 +Counts: 6 +Counts: 6 +Convergence: 0 +Messages: CONVERGENCE: REL_REDUCTION_OF_F <= FACTR*EPSMCH + +CO2 MSE: 4.57 +T MSE: 0.0357 +RMSE: 2.14 + +***Key Metrics*** +TCRE: 1.51 +TCR: 1.77 + +***Historical Warming and ERF*** +GSAT Warming: 0.731 +Ocean Heat Content Change: 471 +Total Aerosol ERF: -1.24 +WMGHG ERF: 3.87 +Methane ERF: 0.54 + +***Future Warming*** +scenario start end GSAT +ssp119 2021 2040 0.73 +ssp119 2041 2060 0.895 +ssp119 2081 2100 0.717 +ssp126 2021 2040 0.745 +ssp126 2041 2060 1.08 +ssp126 2081 2100 1.1 +ssp245 2021 2040 0.75 +ssp245 2041 2060 1.29 +ssp245 2081 2100 1.98 +ssp370 2021 2040 0.761 +ssp370 2041 2060 1.43 +ssp370 2081 2100 2.94 +ssp585 2021 2040 0.877 +ssp585 2041 2060 1.74 +ssp585 2081 2100 3.79 diff --git a/results/03_smooth_temps_5yr.txt b/results/03_smooth_temps_5yr.txt new file mode 100644 index 0000000..ed99482 --- /dev/null +++ b/results/03_smooth_temps_5yr.txt @@ -0,0 +1,43 @@ +parameters values +beta 0.528262971929169 +q10_rh 1.76 +diff 2.37997009562544 + +Objective Function Value: 2.3 +Counts: 6 +Counts: 6 +Convergence: 0 +Messages: CONVERGENCE: REL_REDUCTION_OF_F <= FACTR*EPSMCH + +CO2 MSE: 4.57 +T MSE: 0.0357 +RMSE: 2.14 + +***Key Metrics*** +TCRE: 1.51 +TCR: 1.77 + +***Historical Warming and ERF*** +GSAT Warming: 0.731 +Ocean Heat Content Change: 471 +Total Aerosol ERF: -1.24 +WMGHG ERF: 3.87 +Methane ERF: 0.54 + +***Future Warming*** +scenario start end GSAT +ssp119 2021 2040 0.73 +ssp119 2041 2060 0.895 +ssp119 2081 2100 0.717 +ssp126 2021 2040 0.745 +ssp126 2041 2060 1.08 +ssp126 2081 2100 1.1 +ssp245 2021 2040 0.75 +ssp245 2041 2060 1.29 +ssp245 2081 2100 1.98 +ssp370 2021 2040 0.761 +ssp370 2041 2060 1.43 +ssp370 2081 2100 2.94 +ssp585 2021 2040 0.877 +ssp585 2041 2060 1.74 +ssp585 2081 2100 3.79 diff --git a/results/04_smooth_temps_10yr.txt b/results/04_smooth_temps_10yr.txt new file mode 100644 index 0000000..9ce92b0 --- /dev/null +++ b/results/04_smooth_temps_10yr.txt @@ -0,0 +1,43 @@ +parameters values +beta 0.5282627568606 +q10_rh 1.76 +diff 2.37996967686392 + +Objective Function Value: 2.3 +Counts: 6 +Counts: 6 +Convergence: 0 +Messages: CONVERGENCE: REL_REDUCTION_OF_F <= FACTR*EPSMCH + +CO2 MSE: 4.57 +T MSE: 0.0357 +RMSE: 2.14 + +***Key Metrics*** +TCRE: 1.51 +TCR: 1.77 + +***Historical Warming and ERF*** +GSAT Warming: 0.731 +Ocean Heat Content Change: 471 +Total Aerosol ERF: -1.24 +WMGHG ERF: 3.87 +Methane ERF: 0.54 + +***Future Warming*** +scenario start end GSAT +ssp119 2021 2040 0.73 +ssp119 2041 2060 0.895 +ssp119 2081 2100 0.717 +ssp126 2021 2040 0.745 +ssp126 2041 2060 1.08 +ssp126 2081 2100 1.1 +ssp245 2021 2040 0.75 +ssp245 2041 2060 1.29 +ssp245 2081 2100 1.98 +ssp370 2021 2040 0.761 +ssp370 2041 2060 1.43 +ssp370 2081 2100 2.94 +ssp585 2021 2040 0.877 +ssp585 2041 2060 1.74 +ssp585 2081 2100 3.79 diff --git a/results/05_NMSE.txt b/results/05_NMSE.txt new file mode 100644 index 0000000..40461e5 --- /dev/null +++ b/results/05_NMSE.txt @@ -0,0 +1,43 @@ +parameters values +beta 0.268 +q10_rh 2.64 +diff 2.2 + +Objective Function Value: 0.000565 +Counts: 11 +Counts: 11 +Convergence: 0 +Messages: CONVERGENCE: NORM OF PROJECTED GRADIENT <= PGTOL + +CO2 MSE: 173 +T MSE: 0.0211 +RMSE: 13.2 + +***Key Metrics*** +TCRE: 1.72 +TCR: 1.79 + +***Historical Warming and ERF*** +GSAT Warming: 0.925 +Ocean Heat Content Change: 557 +Total Aerosol ERF: -1.24 +WMGHG ERF: 4.38 +Methane ERF: 0.544 + +***Future Warming*** +scenario start end GSAT +ssp119 2021 2040 0.857 +ssp119 2041 2060 1.14 +ssp119 2081 2100 1.11 +ssp126 2021 2040 0.871 +ssp126 2041 2060 1.32 +ssp126 2081 2100 1.53 +ssp245 2021 2040 0.875 +ssp245 2041 2060 1.53 +ssp245 2081 2100 2.48 +ssp370 2021 2040 0.885 +ssp370 2041 2060 1.67 +ssp370 2081 2100 3.42 +ssp585 2021 2040 1 +ssp585 2041 2060 1.99 +ssp585 2081 2100 4.29 diff --git a/results/05_NMSE_big_box.txt b/results/05_NMSE_big_box.txt new file mode 100644 index 0000000..9962f58 --- /dev/null +++ b/results/05_NMSE_big_box.txt @@ -0,0 +1,43 @@ +parameters values +beta 0 +q10_rh 1.49599440204924 +diff 2.6 + +Objective Function Value: 0.000406 +Counts: 23 +Counts: 23 +Convergence: 0 +Messages: CONVERGENCE: REL_REDUCTION_OF_F <= FACTR*EPSMCH + +CO2 MSE: 918 +T MSE: 0.0143 +RMSE: 30.3 + +***Key Metrics*** +TCRE: 1.69 +TCR: 1.74 + +***Historical Warming and ERF*** +GSAT Warming: 1.1 +Ocean Heat Content Change: 684 +Total Aerosol ERF: -1.24 +WMGHG ERF: 4.88 +Methane ERF: 0.549 + +***Future Warming*** +scenario start end GSAT +ssp119 2021 2040 0.885 +ssp119 2041 2060 1.19 +ssp119 2081 2100 1.17 +ssp126 2021 2040 0.897 +ssp126 2041 2060 1.38 +ssp126 2081 2100 1.61 +ssp245 2021 2040 0.9 +ssp245 2041 2060 1.57 +ssp245 2081 2100 2.53 +ssp370 2021 2040 0.909 +ssp370 2041 2060 1.7 +ssp370 2081 2100 3.42 +ssp585 2021 2040 1.02 +ssp585 2041 2060 2 +ssp585 2081 2100 4.24 diff --git a/results/06_smooth_temps_nmse_3yr.txt b/results/06_smooth_temps_nmse_3yr.txt new file mode 100644 index 0000000..a0997c0 --- /dev/null +++ b/results/06_smooth_temps_nmse_3yr.txt @@ -0,0 +1,43 @@ +parameters values +beta 0.268 +q10_rh 2.64 +diff 2.2 + +Objective Function Value: 0.00046 +Counts: 11 +Counts: 11 +Convergence: 0 +Messages: CONVERGENCE: NORM OF PROJECTED GRADIENT <= PGTOL + +CO2 MSE: 173 +T MSE: 0.0211 +RMSE: 13.2 + +***Key Metrics*** +TCRE: 1.72 +TCR: 1.79 + +***Historical Warming and ERF*** +GSAT Warming: 0.925 +Ocean Heat Content Change: 557 +Total Aerosol ERF: -1.24 +WMGHG ERF: 4.38 +Methane ERF: 0.544 + +***Future Warming*** +scenario start end GSAT +ssp119 2021 2040 0.857 +ssp119 2041 2060 1.14 +ssp119 2081 2100 1.11 +ssp126 2021 2040 0.871 +ssp126 2041 2060 1.32 +ssp126 2081 2100 1.53 +ssp245 2021 2040 0.875 +ssp245 2041 2060 1.53 +ssp245 2081 2100 2.48 +ssp370 2021 2040 0.885 +ssp370 2041 2060 1.67 +ssp370 2081 2100 3.42 +ssp585 2021 2040 1 +ssp585 2041 2060 1.99 +ssp585 2081 2100 4.29 diff --git a/results/06_smooth_temps_nmse_3yr_big_box.txt b/results/06_smooth_temps_nmse_3yr_big_box.txt new file mode 100644 index 0000000..66ded55 --- /dev/null +++ b/results/06_smooth_temps_nmse_3yr_big_box.txt @@ -0,0 +1,43 @@ +parameters values +beta 0 +q10_rh 1.57851439895725 +diff 2.6 + +Objective Function Value: 0.000289 +Counts: 19 +Counts: 19 +Convergence: 0 +Messages: CONVERGENCE: REL_REDUCTION_OF_F <= FACTR*EPSMCH + +CO2 MSE: 937 +T MSE: 0.0143 +RMSE: 30.6 + +***Key Metrics*** +TCRE: 1.7 +TCR: 1.74 + +***Historical Warming and ERF*** +GSAT Warming: 1.1 +Ocean Heat Content Change: 687 +Total Aerosol ERF: -1.24 +WMGHG ERF: 4.89 +Methane ERF: 0.549 + +***Future Warming*** +scenario start end GSAT +ssp119 2021 2040 0.889 +ssp119 2041 2060 1.2 +ssp119 2081 2100 1.2 +ssp126 2021 2040 0.902 +ssp126 2041 2060 1.38 +ssp126 2081 2100 1.63 +ssp245 2021 2040 0.905 +ssp245 2041 2060 1.58 +ssp245 2081 2100 2.56 +ssp370 2021 2040 0.913 +ssp370 2041 2060 1.71 +ssp370 2081 2100 3.44 +ssp585 2021 2040 1.03 +ssp585 2041 2060 2.01 +ssp585 2081 2100 4.26 diff --git a/results/08_smooth_temps_nmse_10yr.txt b/results/08_smooth_temps_nmse_10yr.txt new file mode 100644 index 0000000..0f3332a --- /dev/null +++ b/results/08_smooth_temps_nmse_10yr.txt @@ -0,0 +1,43 @@ +parameters values +beta 0.268 +q10_rh 2.64 +diff 2.2 + +Objective Function Value: 0.000411 +Counts: 11 +Counts: 11 +Convergence: 0 +Messages: CONVERGENCE: NORM OF PROJECTED GRADIENT <= PGTOL + +CO2 MSE: 173 +T MSE: 0.0211 +RMSE: 13.2 + +***Key Metrics*** +TCRE: 1.72 +TCR: 1.79 + +***Historical Warming and ERF*** +GSAT Warming: 0.925 +Ocean Heat Content Change: 557 +Total Aerosol ERF: -1.24 +WMGHG ERF: 4.38 +Methane ERF: 0.544 + +***Future Warming*** +scenario start end GSAT +ssp119 2021 2040 0.857 +ssp119 2041 2060 1.14 +ssp119 2081 2100 1.11 +ssp126 2021 2040 0.871 +ssp126 2041 2060 1.32 +ssp126 2081 2100 1.53 +ssp245 2021 2040 0.875 +ssp245 2041 2060 1.53 +ssp245 2081 2100 2.48 +ssp370 2021 2040 0.885 +ssp370 2041 2060 1.67 +ssp370 2081 2100 3.42 +ssp585 2021 2040 1 +ssp585 2041 2060 1.99 +ssp585 2081 2100 4.29 diff --git a/results/08_smooth_temps_nmse_10yr_big_box.txt b/results/08_smooth_temps_nmse_10yr_big_box.txt new file mode 100644 index 0000000..1a116a4 --- /dev/null +++ b/results/08_smooth_temps_nmse_10yr_big_box.txt @@ -0,0 +1,43 @@ +parameters values +beta 0 +q10_rh 1.95459711198892 +diff 2.6 + +Objective Function Value: 0.000213 +Counts: 19 +Counts: 19 +Convergence: 0 +Messages: CONVERGENCE: REL_REDUCTION_OF_F <= FACTR*EPSMCH + +CO2 MSE: 1020 +T MSE: 0.0143 +RMSE: 31.9 + +***Key Metrics*** +TCRE: 1.73 +TCR: 1.74 + +***Historical Warming and ERF*** +GSAT Warming: 1.12 +Ocean Heat Content Change: 697 +Total Aerosol ERF: -1.24 +WMGHG ERF: 4.95 +Methane ERF: 0.55 + +***Future Warming*** +scenario start end GSAT +ssp119 2021 2040 0.907 +ssp119 2041 2060 1.24 +ssp119 2081 2100 1.29 +ssp126 2021 2040 0.92 +ssp126 2041 2060 1.42 +ssp126 2081 2100 1.72 +ssp245 2021 2040 0.922 +ssp245 2041 2060 1.61 +ssp245 2081 2100 2.64 +ssp370 2021 2040 0.93 +ssp370 2041 2060 1.74 +ssp370 2081 2100 3.51 +ssp585 2021 2040 1.04 +ssp585 2041 2060 2.04 +ssp585 2081 2100 4.33 diff --git a/results/09_unc_nmse.txt b/results/09_unc_nmse.txt new file mode 100644 index 0000000..2e9c578 --- /dev/null +++ b/results/09_unc_nmse.txt @@ -0,0 +1,44 @@ +parameters values +beta 0.268 +q10_rh 2.64 +diff 2.2 + +Objective Function Value: 0.000164 +Counts: 10 +Counts: 10 +Convergence: 0 +Messages: CONVERGENCE: NORM OF PROJECTED GRADIENT <= PGTOL + +CO2 MSE: 173 +T MSE: 0.0211 +RMSE: 13.2 +T MSE accounting for unc: 0.00598 + +***Key Metrics*** +TCRE: 1.72 +TCR: 1.79 + +***Historical Warming and ERF*** +GSAT Warming: 0.925 +Ocean Heat Content Change: 557 +Total Aerosol ERF: -1.24 +WMGHG ERF: 4.38 +Methane ERF: 0.544 + +***Future Warming*** +scenario start end GSAT +ssp119 2021 2040 0.857 +ssp119 2041 2060 1.14 +ssp119 2081 2100 1.11 +ssp126 2021 2040 0.871 +ssp126 2041 2060 1.32 +ssp126 2081 2100 1.53 +ssp245 2021 2040 0.875 +ssp245 2041 2060 1.53 +ssp245 2081 2100 2.48 +ssp370 2021 2040 0.885 +ssp370 2041 2060 1.67 +ssp370 2081 2100 3.42 +ssp585 2021 2040 1 +ssp585 2041 2060 1.99 +ssp585 2081 2100 4.29 diff --git a/results/09_unc_nmse_big_box.txt b/results/09_unc_nmse_big_box.txt new file mode 100644 index 0000000..e654969 --- /dev/null +++ b/results/09_unc_nmse_big_box.txt @@ -0,0 +1,44 @@ +parameters values +beta 0.0275279195575259 +q10_rh 1.75689544007986 +diff 2.6 + +Objective Function Value: 0.000123 +Counts: 12 +Counts: 12 +Convergence: 0 +Messages: CONVERGENCE: REL_REDUCTION_OF_F <= FACTR*EPSMCH + +CO2 MSE: 824 +T MSE: 0.0145 +RMSE: 28.7 +T MSE accounting for unc: 0.00371 + +***Key Metrics*** +TCRE: 1.7 +TCR: 1.74 + +***Historical Warming and ERF*** +GSAT Warming: 1.08 +Ocean Heat Content Change: 676 +Total Aerosol ERF: -1.24 +WMGHG ERF: 4.85 +Methane ERF: 0.549 + +***Future Warming*** +scenario start end GSAT +ssp119 2021 2040 0.888 +ssp119 2041 2060 1.2 +ssp119 2081 2100 1.21 +ssp126 2021 2040 0.9 +ssp126 2041 2060 1.38 +ssp126 2081 2100 1.64 +ssp245 2021 2040 0.903 +ssp245 2041 2060 1.58 +ssp245 2081 2100 2.56 +ssp370 2021 2040 0.912 +ssp370 2041 2060 1.71 +ssp370 2081 2100 3.44 +ssp585 2021 2040 1.02 +ssp585 2041 2060 2.01 +ssp585 2081 2100 4.27 diff --git a/results/10_ecs_unc_nmse.txt b/results/10_ecs_unc_nmse.txt new file mode 100644 index 0000000..96079f8 --- /dev/null +++ b/results/10_ecs_unc_nmse.txt @@ -0,0 +1,45 @@ +parameters values +beta 0.268 +q10_rh 2.64 +diff 2.4 +S 3.97216492744402 + +Objective Function Value: 0.000138 +Counts: 19 +Counts: 19 +Convergence: 0 +Messages: CONVERGENCE: REL_REDUCTION_OF_F <= FACTR*EPSMCH + +CO2 MSE: 182 +T MSE: 0.017 +RMSE: 13.5 +T MSE accounting for unc: 0.00501 + +***Key Metrics*** +TCRE: 2.16 +TCR: 2.07 + +***Historical Warming and ERF*** +GSAT Warming: 1.08 +Ocean Heat Content Change: 659 +Total Aerosol ERF: -1.24 +WMGHG ERF: 4.41 +Methane ERF: 0.548 + +***Future Warming*** +scenario start end GSAT +ssp119 2021 2040 1.02 +ssp119 2041 2060 1.42 +ssp119 2081 2100 1.55 +ssp126 2021 2040 1.03 +ssp126 2041 2060 1.62 +ssp126 2081 2100 2.04 +ssp245 2021 2040 1.04 +ssp245 2041 2060 1.84 +ssp245 2081 2100 3.11 +ssp370 2021 2040 1.05 +ssp370 2041 2060 1.99 +ssp370 2081 2100 4.15 +ssp585 2021 2040 1.17 +ssp585 2041 2060 2.35 +ssp585 2081 2100 5.15 diff --git a/results/10_ecs_unc_nmse_big_box.txt b/results/10_ecs_unc_nmse_big_box.txt new file mode 100644 index 0000000..7af9b27 --- /dev/null +++ b/results/10_ecs_unc_nmse_big_box.txt @@ -0,0 +1,45 @@ +parameters values +beta 0.00582111555564215 +q10_rh 1.00065113643678 +diff 2.6 +S 3.16051323043557 + +Objective Function Value: 0.000122 +Counts: 43 +Counts: 43 +Convergence: 0 +Messages: CONVERGENCE: REL_REDUCTION_OF_F <= FACTR*EPSMCH + +CO2 MSE: 762 +T MSE: 0.0145 +RMSE: 27.6 +T MSE accounting for unc: 0.00372 + +***Key Metrics*** +TCRE: 1.68 +TCR: 1.79 + +***Historical Warming and ERF*** +GSAT Warming: 1.09 +Ocean Heat Content Change: 682 +Total Aerosol ERF: -1.24 +WMGHG ERF: 4.76 +Methane ERF: 0.549 + +***Future Warming*** +scenario start end GSAT +ssp119 2021 2040 0.88 +ssp119 2041 2060 1.17 +ssp119 2081 2100 1.05 +ssp126 2021 2040 0.893 +ssp126 2041 2060 1.35 +ssp126 2081 2100 1.51 +ssp245 2021 2040 0.897 +ssp245 2041 2060 1.55 +ssp245 2081 2100 2.48 +ssp370 2021 2040 0.906 +ssp370 2041 2060 1.69 +ssp370 2081 2100 3.4 +ssp585 2021 2040 1.02 +ssp585 2041 2060 2 +ssp585 2081 2100 4.26 diff --git a/results/11_alpha_ecs_unc_nmse.txt b/results/11_alpha_ecs_unc_nmse.txt new file mode 100644 index 0000000..3067834 --- /dev/null +++ b/results/11_alpha_ecs_unc_nmse.txt @@ -0,0 +1,46 @@ +parameters values +beta 0.569553577866329 +q10_rh 1.76 +diff 2.38453856930776 +S 2.95650938955361 +alpha 0.49168059960629 + +Objective Function Value: 6.75e-05 +Counts: 11 +Counts: 11 +Convergence: 0 +Messages: CONVERGENCE: REL_REDUCTION_OF_F <= FACTR*EPSMCH + +CO2 MSE: 5.03 +T MSE: 0.0129 +RMSE: 2.24 +T MSE accounting for unc: 0.00254 + +***Key Metrics*** +TCRE: 1.48 +TCR: 1.75 + +***Historical Warming and ERF*** +GSAT Warming: 1.07 +Ocean Heat Content Change: 623 +Total Aerosol ERF: -1.24 +WMGHG ERF: 3.88 +Methane ERF: 0.547 + +***Future Warming*** +scenario start end GSAT +ssp119 2021 2040 0.611 +ssp119 2041 2060 0.683 +ssp119 2081 2100 0.426 +ssp126 2021 2040 0.642 +ssp126 2041 2060 0.881 +ssp126 2081 2100 0.8 +ssp245 2021 2040 0.692 +ssp245 2041 2060 1.18 +ssp245 2081 2100 1.76 +ssp370 2021 2040 0.735 +ssp370 2041 2060 1.4 +ssp370 2081 2100 2.87 +ssp585 2021 2040 0.786 +ssp585 2041 2060 1.59 +ssp585 2081 2100 3.57 diff --git a/results/11_alpha_ecs_unc_nmse_big_box.txt b/results/11_alpha_ecs_unc_nmse_big_box.txt new file mode 100644 index 0000000..34cc2cb --- /dev/null +++ b/results/11_alpha_ecs_unc_nmse_big_box.txt @@ -0,0 +1,46 @@ +parameters values +beta 0.50229087166652 +q10_rh 0.986962103838984 +diff 2 +S 2.88436966557328 +alpha 0.499997310177854 + +Objective Function Value: 6.62e-05 +Counts: 40 +Counts: 40 +Convergence: 0 +Messages: CONVERGENCE: REL_REDUCTION_OF_F <= FACTR*EPSMCH + +CO2 MSE: 4.86 +T MSE: 0.0127 +RMSE: 2.2 +T MSE accounting for unc: 0.00249 + +***Key Metrics*** +TCRE: 1.33 +TCR: 1.78 + +***Historical Warming and ERF*** +GSAT Warming: 1.07 +Ocean Heat Content Change: 582 +Total Aerosol ERF: -1.24 +WMGHG ERF: 3.84 +Methane ERF: 0.547 + +***Future Warming*** +scenario start end GSAT +ssp119 2021 2040 0.587 +ssp119 2041 2060 0.616 +ssp119 2081 2100 0.26 +ssp126 2021 2040 0.618 +ssp126 2041 2060 0.821 +ssp126 2081 2100 0.638 +ssp245 2021 2040 0.669 +ssp245 2041 2060 1.12 +ssp245 2081 2100 1.61 +ssp370 2021 2040 0.714 +ssp370 2041 2060 1.35 +ssp370 2081 2100 2.74 +ssp585 2021 2040 0.767 +ssp585 2041 2060 1.55 +ssp585 2081 2100 3.45 diff --git a/results/12_alpha_ecs_ohc_unc_nmse.txt b/results/12_alpha_ecs_ohc_unc_nmse.txt new file mode 100644 index 0000000..62f6bac --- /dev/null +++ b/results/12_alpha_ecs_ohc_unc_nmse.txt @@ -0,0 +1,47 @@ +parameters values +beta 0.64953154546281 +q10_rh 1.76 +diff 1.042 +S 2.32537207553702 +alpha 0.438414955305681 + +Objective Function Value: 4.77e-05 +Counts: 40 +Counts: 40 +Convergence: 0 +Messages: CONVERGENCE: REL_REDUCTION_OF_F <= FACTR*EPSMCH + +CO2 MSE: 8.57 +T MSE: 0.0127 +RMSE: 2.93 +T MSE accounting for unc: 0.0026 +OHC MSE accounting for unc: 9.71 + +***Key Metrics*** +TCRE: 1.23 +TCR: 1.67 + +***Historical Warming and ERF*** +GSAT Warming: 1.01 +Ocean Heat Content Change: 425 +Total Aerosol ERF: -1.24 +WMGHG ERF: 3.78 +Methane ERF: 0.545 + +***Future Warming*** +scenario start end GSAT +ssp119 2021 2040 0.51 +ssp119 2041 2060 0.484 +ssp119 2081 2100 0.138 +ssp126 2021 2040 0.545 +ssp126 2041 2060 0.688 +ssp126 2081 2100 0.473 +ssp245 2021 2040 0.603 +ssp245 2041 2060 1 +ssp245 2081 2100 1.4 +ssp370 2021 2040 0.654 +ssp370 2041 2060 1.24 +ssp370 2081 2100 2.51 +ssp585 2021 2040 0.7 +ssp585 2041 2060 1.43 +ssp585 2081 2100 3.16 diff --git a/results/13_alpha_ecs_ohc_mvsse.txt b/results/13_alpha_ecs_ohc_mvsse.txt new file mode 100644 index 0000000..0e1514c --- /dev/null +++ b/results/13_alpha_ecs_ohc_mvsse.txt @@ -0,0 +1,48 @@ +parameters values +beta 0.526231809997589 +q10_rh 2.31055396561195 +diff 1.042 +S 2.82788331898928 +alpha 1.40507672116753 + +Objective Function Value: 91.3 +Counts: 104 +Counts: 104 +Convergence: 0 +Messages: CONVERGENCE: REL_REDUCTION_OF_F <= FACTR*EPSMCH + +CO2 MSE: 3.4 +T MSE: 0.102 +RMSE: 1.84 +CO2 MVSSE: 236 +T MVSSE: 37.3 +OHC MVSSE: 0.385 + +***Key Metrics*** +TCRE: 1.62 +TCR: 1.92 + +***Historical Warming and ERF*** +GSAT Warming: 0.502 +Ocean Heat Content Change: 291 +Total Aerosol ERF: -1.24 +WMGHG ERF: 3.88 +Methane ERF: 0.536 + +***Future Warming*** +scenario start end GSAT +ssp119 2021 2040 0.904 +ssp119 2041 2060 1.13 +ssp119 2081 2100 0.968 +ssp126 2021 2040 0.907 +ssp126 2041 2060 1.33 +ssp126 2081 2100 1.36 +ssp245 2021 2040 0.872 +ssp245 2041 2060 1.49 +ssp245 2081 2100 2.3 +ssp370 2021 2040 0.854 +ssp370 2041 2060 1.57 +ssp370 2081 2100 3.23 +ssp585 2021 2040 1.05 +ssp585 2041 2060 2.04 +ssp585 2081 2100 4.3 diff --git a/results/alpha_comparison_plots.jpeg b/results/alpha_comparison_plots.jpeg new file mode 100644 index 0000000..28522c4 Binary files /dev/null and b/results/alpha_comparison_plots.jpeg differ diff --git a/results/default_stats.txt b/results/default_stats.txt new file mode 100644 index 0000000..7f93983 --- /dev/null +++ b/results/default_stats.txt @@ -0,0 +1,60 @@ +scenario year variable value units +Unnamed Hector core NA beta 0.53 (unitless) +Unnamed Hector core NA q10_rh 1.76 (unitless) +Unnamed Hector core NA diff 2.38 cm2/s + +Mean of T, CO2 MSEs: 2.31 +Mean of T, CO2 NMSEs: 0.000948 +Mean of T, CO2 NMSEs (with unc): 0.000313 +Mean of T, CO2, OHC NMSEs (with unc): 0.00022 +Mean of T, CO2, OHC MVSSEs: 113 + +Mean of smooth T, CO2 MSEs (k = 3 ): 2.3 +Mean of smooth T, CO2 NMSEs (k = 3 ): 0.000864 +Mean of smooth T, CO2 MSEs (k = 5 ): 2.3 +Mean of smooth T, CO2 NMSEs (k = 5 ): 0.000846 +Mean of smooth T, CO2 MSEs (k = 10 ): 2.3 +Mean of smooth T, CO2 NMSEs (k = 10 ): 0.000836 + +CO2 MSE: 4.58 +T MSE: 0.0358 +T MSE with unc: 0.0118 +RMSE: 2.14 + +CO2 NMSE: 2.74e-07 +T NMSE: 0.0019 +T NMSE with unc: 0.000625 +OHC NMSE with unc: 3.43e-05 + +CO2 MVSSE: 318 +T MVSSE: 19.5 +OHC MVSSE: 1.02 + +***Key Metrics*** +TCRE: 1.51 +TCR: 1.77 + +***Historical Warming and ERF*** +GSAT Warming: 0.73 +Ocean Heat Content Change: 471 +Total Aerosol ERF: -1.24 +WMGHG ERF: 3.87 +Methane ERF: 0.54 + +***Future Warming*** +scenario start end GSAT +ssp119 2021 2040 0.73 +ssp119 2041 2060 0.894 +ssp119 2081 2100 0.716 +ssp126 2021 2040 0.744 +ssp126 2041 2060 1.08 +ssp126 2081 2100 1.1 +ssp245 2021 2040 0.75 +ssp245 2041 2060 1.29 +ssp245 2081 2100 1.98 +ssp370 2021 2040 0.761 +ssp370 2041 2060 1.43 +ssp370 2081 2100 2.94 +ssp585 2021 2040 0.877 +ssp585 2041 2060 1.74 +ssp585 2081 2100 3.79 diff --git a/results/nmse_comparison_plots.jpeg b/results/nmse_comparison_plots.jpeg new file mode 100644 index 0000000..950198a Binary files /dev/null and b/results/nmse_comparison_plots.jpeg differ diff --git a/results/nmse_unc_comparison_plots.jpeg b/results/nmse_unc_comparison_plots.jpeg new file mode 100644 index 0000000..688b5b5 Binary files /dev/null and b/results/nmse_unc_comparison_plots.jpeg differ diff --git a/results/ohc_plot.jpeg b/results/ohc_plot.jpeg new file mode 100644 index 0000000..e685898 Binary files /dev/null and b/results/ohc_plot.jpeg differ diff --git a/scripts/alpha_ohc_mvsse.R b/scripts/alpha_ohc_mvsse.R new file mode 100644 index 0000000..4569ae2 --- /dev/null +++ b/scripts/alpha_ohc_mvsse.R @@ -0,0 +1,96 @@ +# Script for experiment 13 +# Script to optimize over T, OHC and CO2 MVSSEs +# Also includes ECS and alpha as params to optimize over +# Uses OHC range from Matilda table rather than manuscript script +# Author: Peter Scully +# Date: 6/25/24 + +### Constants and Imports ### + +# Importing libraries +library(hector) + +# Setting up file paths +COMP_DATA_DIR <- file.path(here::here(), "comparison_data") +SCRIPTS_DIR <- file.path(here::here(), "scripts") +RESULTS_DIR <- file.path(here::here(), "results") + +CO2_PATH <- file.path(COMP_DATA_DIR, + "Supplementary_Table_UoM_GHGConcentrations-1-1-0_annualmeans_v23March2017.csv") +TEMP_PATH <- + file.path(COMP_DATA_DIR, + "HadCRUT.5.0.2.0.analysis.summary_series.global.annual.csv") + +OHC_PATH <- file.path(COMP_DATA_DIR, "OHC_ensemble_Kuhlbrodt_etal_2022.csv") + + +INI_FILE <- system.file("input/hector_ssp245.ini", package = "hector") +PARAMS <- c(BETA(), Q10_RH(), DIFFUSIVITY(), ECS(), AERO_SCALE()) + +OUTPUT <- file.path(RESULTS_DIR, "13_alpha_ecs_ohc_mvsse.txt") + + +source(file.path(SCRIPTS_DIR, "major_functions.R")) + +### Getting observational data ### +co2_data <- get_co2_data(CO2_PATH, include_unc = TRUE) +temp_data <- get_temp_data(TEMP_PATH, include_unc = TRUE) +ohc_data <- get_ohc_data(OHC_PATH, include_unc = T) +obs_data <- rbind(co2_data, temp_data, ohc_data) + +### Calling optim ### +best_pars <- run_optim(obs_data = obs_data, + ini_file = INI_FILE, + params = PARAMS, + lower = c(0.5 - 0.232, 2.2 - 0.44, 1.16 - 0.118, 2, 0), + upper = c(0.5 + 0.232, 2.2 + 0.44, 1.16 + 0.118, 5, 3), + yrs = 1750:2014, + vars = c(GMST(), CONCENTRATIONS_CO2(), HEAT_FLUX()), + error_fn = mean_T_CO2_OHC_mvsse, + include_unc = T, + method = "L-BFGS-B", + output_file = OUTPUT) + +### Outputting individual MSEs ### +hector_data <- run_hector(ini_file = INI_FILE, + params = PARAMS, + vals = best_pars, + yrs = 1750:2014, + vars = c(GMST(), CONCENTRATIONS_CO2(), HEAT_FLUX())) + +T_mse <- get_var_mse(obs_data = obs_data, + hector_data = hector_data, + var = GMST(), + yrs = 1850:2014) +CO2_mse <- get_var_mse(obs_data = obs_data, + hector_data = hector_data, + var = CONCENTRATIONS_CO2(), + yrs = c(1750, 1850:2014)) + +# Getting MVSSEs +T_mvsse <- get_var_mvsse(obs_data = obs_data, + hector_data = hector_data, + var = GMST(), + yrs = 1850:2014, + mse_fn = mvsse) +CO2_mvsse <- get_var_mvsse(obs_data = obs_data, + hector_data = hector_data, + var = CONCENTRATIONS_CO2(), + yrs = c(1750, 1850:2014), + mse_fn = mvsse) +OHC_mvsse <- get_var_mvsse(obs_data = obs_data, + hector_data = hector_data, + var = "OHC", + yrs = 1957:2014, + mse_fn = mvsse) + +write_metric("CO2 MSE:", CO2_mse, OUTPUT) +write_metric("T MSE: ", T_mse, OUTPUT) +write_metric("RMSE: ", sqrt(mean(CO2_mse, T_mse)), OUTPUT) # not 100% sure this is how we want to calculate this +write_metric("CO2 MVSSE:", CO2_mvsse, OUTPUT) +write_metric("T MVSSE: ", T_mvsse, OUTPUT) +write_metric("OHC MVSSE:", OHC_mvsse, OUTPUT) +write("", OUTPUT, append = TRUE) + +### Outputting table metrics ### +calc_table_metrics(PARAMS, best_pars, OUTPUT) \ No newline at end of file diff --git a/scripts/alpha_ohc_unc_nmse.R b/scripts/alpha_ohc_unc_nmse.R new file mode 100644 index 0000000..1eab2cf --- /dev/null +++ b/scripts/alpha_ohc_unc_nmse.R @@ -0,0 +1,88 @@ +# Script for experiment 12 +# Script to use normalized T, OHC and CO2 MSEs while accounting for T, OHC unc +# Also includes ECS and alpha as params to optimize over +# Uses OHC range from Matilda table rather than manuscript script +# Author: Peter Scully +# Date: 6/25/24 + +### Constants and Imports ### + +# Importing libraries +library(hector) + +# Setting up file paths +COMP_DATA_DIR <- file.path(here::here(), "comparison_data") +SCRIPTS_DIR <- file.path(here::here(), "scripts") +RESULTS_DIR <- file.path(here::here(), "results") + +CO2_PATH <- file.path(COMP_DATA_DIR, + "Supplementary_Table_UoM_GHGConcentrations-1-1-0_annualmeans_v23March2017.csv") +TEMP_PATH <- + file.path(COMP_DATA_DIR, + "HadCRUT.5.0.2.0.analysis.summary_series.global.annual.csv") + +OHC_PATH <- file.path(COMP_DATA_DIR, "OHC_ensemble_Kuhlbrodt_etal_2022.csv") + + +INI_FILE <- system.file("input/hector_ssp245.ini", package = "hector") +PARAMS <- c(BETA(), Q10_RH(), DIFFUSIVITY(), ECS(), AERO_SCALE()) + +OUTPUT <- file.path(RESULTS_DIR, "12_alpha_ecs_ohc_unc_nmse.txt") + + +source(file.path(SCRIPTS_DIR, "major_functions.R")) + +### Getting observational data ### +co2_data <- get_co2_data(CO2_PATH, include_unc = TRUE) +temp_data <- get_temp_data(TEMP_PATH, include_unc = TRUE) +ohc_data <- get_ohc_data(OHC_PATH, include_unc = T) +obs_data <- rbind(co2_data, temp_data, ohc_data) + +### Calling optim ### +best_pars <- run_optim(obs_data = obs_data, + ini_file = INI_FILE, + params = PARAMS, + lower = c(0.5 - 0.232, 2.2 - 0.44, 1.16 - 0.118, 2, 0), + upper = c(0.5 + 0.232, 2.2 + 0.44, 1.16 + 0.118, 5, 3), + yrs = 1750:2014, + vars = c(GMST(), CONCENTRATIONS_CO2(), HEAT_FLUX()), + error_fn = mean_T_CO2_OHC_nmse_unc, + include_unc = T, + method = "L-BFGS-B", + output_file = OUTPUT) + +### Outputting individual MSEs ### +hector_data <- run_hector(ini_file = INI_FILE, + params = PARAMS, + vals = best_pars, + yrs = 1750:2014, + vars = c(GMST(), CONCENTRATIONS_CO2(), HEAT_FLUX())) + +T_mse <- get_var_mse(obs_data = obs_data, + hector_data = hector_data, + var = GMST(), + yrs = 1850:2014) +CO2_mse <- get_var_mse(obs_data = obs_data, + hector_data = hector_data, + var = CONCENTRATIONS_CO2(), + yrs = c(1750, 1850:2014)) +T_mse_unc <- get_var_mse_unc(obs_data = obs_data, + hector_data = hector_data, + var = GMST(), + yrs = 1850:2014, + mse_fn = mse_unc) +OHC_mse_unc <- get_var_mse_unc(obs_data = obs_data, + hector_data = hector_data, + var = "OHC", + yrs = 1957:2014, + mse_fn = mse_unc) + +write_metric("CO2 MSE:", CO2_mse, OUTPUT) +write_metric("T MSE: ", T_mse, OUTPUT) +write_metric("RMSE: ", sqrt(mean(CO2_mse, T_mse)), OUTPUT) # not 100% sure this is how we want to calculate this +write_metric("T MSE accounting for unc:", T_mse_unc, OUTPUT) +write_metric("OHC MSE accounting for unc:", OHC_mse_unc, OUTPUT) +write("", OUTPUT, append = TRUE) + +### Outputting table metrics ### +calc_table_metrics(PARAMS, best_pars, OUTPUT) \ No newline at end of file diff --git a/scripts/default_stats.R b/scripts/default_stats.R index e5f2473..46e10d8 100644 --- a/scripts/default_stats.R +++ b/scripts/default_stats.R @@ -19,6 +19,8 @@ TEMP_PATH <- file.path(COMP_DATA_DIR, "HadCRUT.5.0.2.0.analysis.summary_series.global.annual.csv") +OHC_PATH <- file.path(COMP_DATA_DIR, "OHC_ensemble_Kuhlbrodt_etal_2022.csv") + INI_FILE <- system.file("input/hector_ssp245.ini", package = "hector") PARAMS <- c(BETA(), Q10_RH(), DIFFUSIVITY()) @@ -31,7 +33,8 @@ source(file.path(SCRIPTS_DIR, "major_functions.R")) ### Getting observational data ### co2_data <- get_co2_data(CO2_PATH, include_unc = T) temp_data <- get_temp_data(TEMP_PATH, include_unc = T) -obs_data <- rbind(co2_data, temp_data) +ohc_data <- get_ohc_data(OHC_PATH, include_unc = T) +obs_data <- rbind(co2_data, temp_data, ohc_data) ### Outputting default values ### @@ -53,7 +56,7 @@ hector_data <- run_hector(ini_file = INI_FILE, params = NULL, vals = best_pars, yrs = 1750:2014, - vars = c(GMST(), CONCENTRATIONS_CO2()), + vars = c(GMST(), CONCENTRATIONS_CO2(), HEAT_FLUX()), include_unc = T) @@ -67,6 +70,12 @@ write_metric("Mean of T, CO2 NMSEs:", write_metric("Mean of T, CO2 NMSEs (with unc):", mean_T_CO2_nmse_unc(obs_data, hector_data), OUTPUT) +write_metric("Mean of T, CO2, OHC NMSEs (with unc):", + mean_T_CO2_OHC_nmse_unc(obs_data, hector_data), + OUTPUT) +write_metric("Mean of T, CO2, OHC MVSSEs:", + mean_T_CO2_OHC_mvsse(obs_data, hector_data), + OUTPUT) write("", OUTPUT, append = TRUE) # Outputting smoothed MSEs @@ -99,7 +108,7 @@ T_mse <- get_var_mse(obs_data = obs_data, T_mse_unc <- get_var_mse_unc(obs_data = obs_data, hector_data = hector_data, var = GMST(), - yrs = c(1850:2014), + yrs = 1850:2014, mse_fn = mse_unc) CO2_mse <- get_var_mse(obs_data = obs_data, hector_data = hector_data, @@ -115,13 +124,38 @@ T_nmse <- get_var_mse(obs_data = obs_data, T_nmse_unc <- get_var_mse_unc(obs_data = obs_data, hector_data = hector_data, var = GMST(), - yrs = c(1850:2014), + yrs = 1850:2014, mse_fn = nmse_unc) CO2_nmse <- get_var_mse(obs_data = obs_data, hector_data = hector_data, var = CONCENTRATIONS_CO2(), yrs = c(1750, 1850:2014), mse_fn = nmse) +OHC_nmse_unc <- get_var_mse_unc(obs_data = obs_data, + hector_data = hector_data, + var = "OHC", + yrs = 1957:2014, + mse_fn = nmse_unc) + + + + +# Getting MVSSEs +T_mvsse <- get_var_mvsse(obs_data = obs_data, + hector_data = hector_data, + var = GMST(), + yrs = 1850:2014, + mse_fn = mvsse) +CO2_mvsse <- get_var_mvsse(obs_data = obs_data, + hector_data = hector_data, + var = CONCENTRATIONS_CO2(), + yrs = c(1750, 1850:2014), + mse_fn = mvsse) +OHC_mvsse <- get_var_mvsse(obs_data = obs_data, + hector_data = hector_data, + var = "OHC", + yrs = 1957:2014, + mse_fn = mvsse) write_metric("CO2 MSE: ", CO2_mse, OUTPUT) write_metric("T MSE: ", T_mse, OUTPUT) @@ -131,6 +165,11 @@ write("", OUTPUT, append = TRUE) write_metric("CO2 NMSE:", CO2_nmse, OUTPUT) write_metric("T NMSE: ", T_nmse, OUTPUT) write_metric("T NMSE with unc:", T_nmse_unc, OUTPUT) +write_metric("OHC NMSE with unc:", OHC_nmse_unc, OUTPUT) +write("", OUTPUT, append = TRUE) +write_metric("CO2 MVSSE:", CO2_mvsse, OUTPUT) +write_metric("T MVSSE: ", T_mvsse, OUTPUT) +write_metric("OHC MVSSE:", OHC_mvsse, OUTPUT) write("", OUTPUT, append = TRUE) ### Outputting table metrics ### diff --git a/scripts/ecs_unc_nmse.R b/scripts/ecs_unc_nmse.R index 1bf31ce..a9b0440 100644 --- a/scripts/ecs_unc_nmse.R +++ b/scripts/ecs_unc_nmse.R @@ -1,3 +1,4 @@ +# Script for running experiment 10 # Script to use normalized the T and CO2 MSEs while accounting for T uncertainty # Also includes ECS as a param to optimize over # Author: Peter Scully @@ -22,7 +23,7 @@ TEMP_PATH <- INI_FILE <- system.file("input/hector_ssp245.ini", package = "hector") PARAMS <- c(BETA(), Q10_RH(), DIFFUSIVITY(), ECS()) -OUTPUT <- file.path(RESULTS_DIR, "ecs_unc_nmse_big_box.txt") +OUTPUT <- file.path(RESULTS_DIR, "10_ecs_unc_nmse_big_box.txt") source(file.path(SCRIPTS_DIR, "major_functions.R")) diff --git a/scripts/error_functions.R b/scripts/error_functions.R index 7000a9b..1b742cf 100644 --- a/scripts/error_functions.R +++ b/scripts/error_functions.R @@ -188,11 +188,18 @@ get_var_mse_unc <- function(obs_data, hector_data, var, yrs, mse_fn) { # # Note: Assumes observed data contains symmetric upper and lower bounds 1 SD # away from actual value -get_var_mse_unc <- function(obs_data, hector_data, var, yrs, mse_fn) { +get_var_mvsse <- function(obs_data, hector_data, var, yrs, mse_fn) { + + # Getting x and sd x <- filter(obs_data, year %in% yrs & variable == var)$value x_upper <- filter(obs_data, year %in% yrs & variable == var)$upper - sd <- x_upper - x - y <- filter(hector_data, year %in% yrs & variable == var)$value + if (var == GMST()) { + sd <- (x_upper - x) / 1.96 # Accounting for 95% confidence interval + } else { + sd <- x_upper - x + } + + y <- filter(hector_data, year %in% yrs & variable == var)$value return(mvsse(x = x, sd = sd, y = y)) } @@ -269,6 +276,61 @@ mean_T_CO2_nmse_unc <- function(obs_data, hector_data) { return(mean(c(T_mse, CO2_mse))) } +# mean_T_CO2_OHC_nmse_unc: function to find the mean of the temperature, CO2, +# and OHC NMSEs between observed and predicted data +# while accounting for temperature and OHC uncertainty +# +# args: +# obs_data - data frame of observed data formatted like Hector data frame +# hector_data - data frame outputted by Hector +# +# Returns: Average NMSE between predicted and observed data +mean_T_CO2_OHC_nmse_unc <- function(obs_data, hector_data) { + T_mse <- get_var_mse_unc(obs_data = obs_data, + hector_data = hector_data, + var = GMST(), + yrs = 1850:2014, + mse_fn = nmse_unc) + CO2_mse <- get_var_mse(obs_data = obs_data, + hector_data = hector_data, + var = CONCENTRATIONS_CO2(), + yrs = c(1750, 1850:2014), + mse_fn = nmse) + OHC_mse <- get_var_mse_unc(obs_data = obs_data, + hector_data = hector_data, + var = "OHC", + yrs = 1957:2014, + mse_fn = nmse_unc) + return(mean(c(T_mse, CO2_mse, OHC_mse))) +} + +# mean_T_CO2_OHC_mvsse: function to find the mean of the temperature, CO2, +# and OHC MVSSEs between observed and predicted data +# +# args: +# obs_data - data frame of observed data formatted like Hector data frame +# hector_data - data frame outputted by Hector +# +# Returns: Average NMSE between predicted and observed data +mean_T_CO2_OHC_mvsse <- function(obs_data, hector_data) { + T_mse <- get_var_mvsse(obs_data = obs_data, + hector_data = hector_data, + var = GMST(), + yrs = 1850:2014, + mse_fn = mvsse) + CO2_mse <- get_var_mvsse(obs_data = obs_data, + hector_data = hector_data, + var = CONCENTRATIONS_CO2(), + yrs = c(1750, 1850:2014), + mse_fn = mvsse) + OHC_mse <- get_var_mvsse(obs_data = obs_data, + hector_data = hector_data, + var = "OHC", + yrs = 1957:2014, + mse_fn = mvsse) + return(mean(c(T_mse, CO2_mse, OHC_mse))) +} + # smooth_T_CO2_mse: function to find the mean of smoothed temperature and CO2 # MSEs between observed & predicted data for a given variable # diff --git a/scripts/graph_comparison_plots.R b/scripts/graph_comparison_plots.R index 78dfbcd..fbf68be 100644 --- a/scripts/graph_comparison_plots.R +++ b/scripts/graph_comparison_plots.R @@ -22,7 +22,7 @@ TEMP_PATH <- INI_FILE <- system.file("input/hector_ssp245.ini", package = "hector") PARAMS <- c(BETA(), Q10_RH(), DIFFUSIVITY(), ECS(), AERO_SCALE()) -OUTPUT <- file.path(RESULTS_DIR, "alpha_initial_comparison_plots.jpeg") +OUTPUT <- file.path(RESULTS_DIR, "alpha_comparison_plots.jpeg") source(file.path(SCRIPTS_DIR, "major_functions.R")) @@ -46,56 +46,57 @@ default_data <- run_hector(ini_file = INI_FILE, default_data$scenario <- "Hector - Default Fit" -# nmse_data <- run_hector(ini_file = INI_FILE, -# params = PARAMS, -# vals = c(0.732, 2.64, 2.4, 3, 1), -# yrs = 1750:2014, -# vars = c(GMST(), CONCENTRATIONS_CO2())) -# nmse_data$scenario <- "Hector - Fit to NMSEs w/ unc" -# -# nmse_bb_data <- run_hector(ini_file = INI_FILE, -# params = PARAMS, -# vals = c(1.084, 3.52, 2, 3, 1), -# yrs = 1750:2014, -# vars = c(GMST(), CONCENTRATIONS_CO2())) -# nmse_bb_data$scenario <- "Hector - Fit to NMSEs w/ unc, big box" - -nmse_ecs_data <- run_hector(ini_file = INI_FILE, +nmse_data <- run_hector(ini_file = INI_FILE, params = PARAMS, - vals = c(0.732, 2.24, 2.4, 5, 1), + vals = c(0.268, 2.64, 2.2, 3, 1), yrs = 1750:2014, vars = c(GMST(), CONCENTRATIONS_CO2())) -nmse_ecs_data$scenario <- "Hector - NMSEs w/ unc & Tuning S" +nmse_data$scenario <- "Hector - NMSE w/ unc Fit" -nmse_bb_ecs_data <- run_hector(ini_file = INI_FILE, +nmse_bb_data <- run_hector(ini_file = INI_FILE, + params = PARAMS, + vals = c(0.028, 1.76, 2.6, 3, 1), + yrs = 1750:2014, + vars = c(GMST(), CONCENTRATIONS_CO2())) +nmse_bb_data$scenario <- "Hector - NMSE w/ unc, Big Box Fit" + +ecs_data <- run_hector(ini_file = INI_FILE, params = PARAMS, - vals = c(1.069, 3.52, 2, 5, 1), + vals = c(0.268, 1.95, 2.6, 3.97, 1), yrs = 1750:2014, vars = c(GMST(), CONCENTRATIONS_CO2())) -nmse_bb_ecs_data$scenario <- "Hector - NMSEs w/ unc, big box & Tuning S" +ecs_data$scenario <- "Hector - NMSE w/ unc, Fit w/ S" + +ecs_bb_data <- run_hector(ini_file = INI_FILE, + params = PARAMS, + vals = c(0.006, 1, 2.6, 3.16, 1), + yrs = 1750:2014, + vars = c(GMST(), CONCENTRATIONS_CO2())) +ecs_bb_data$scenario <- "Hector - NMSE w/ unc, Big Box Fit w/ S" alpha_data <- run_hector(ini_file = INI_FILE, - params = PARAMS, - vals = c(0.732, 2.24, 2.4, 5, 1.15), - yrs = 1750:2014, - vars = c(GMST(), CONCENTRATIONS_CO2())) -alpha_data$scenario <- "Hector - NMSEs w/ unc & Tuning S and Alpha" + params = PARAMS, + vals = c(0.57, 1.76, 2.38, 2.96, 0.492), + yrs = 1750:2014, + vars = c(GMST(), CONCENTRATIONS_CO2())) +alpha_data$scenario <- "Hector - NMSE w/ unc, Fit w/ S, Alpha" alpha_bb_data <- run_hector(ini_file = INI_FILE, - params = PARAMS, - vals = c(1.196, 3.52, 2, 5, 0.948), - yrs = 1750:2014, - vars = c(GMST(), CONCENTRATIONS_CO2())) -alpha_bb_data$scenario <- "Hector - NMSEs w/ unc, big box & Tuning S and Alpha" - -alpha_vbb_data <- run_hector(ini_file = INI_FILE, params = PARAMS, - vals = c(1.66, 4.4, 1.8, 6, 0.83), + vals = c(0.502, 0.99, 2, 2.88, 0.5), yrs = 1750:2014, vars = c(GMST(), CONCENTRATIONS_CO2())) -alpha_vbb_data$scenario <- "Hector - NMSEs w/ unc, very big box & Tuning S and Alpha" +alpha_bb_data$scenario <- "Hector - NMSE w/ unc, Big Box Fit w/ S, Alpha" + +alpha_ohc_data <- run_hector(ini_file = INI_FILE, + params = PARAMS, + vals = c(0.65, 1.76, 1.04, 2.33, 0.438), + yrs = 1750:2014, + vars = c(GMST(), CONCENTRATIONS_CO2())) +alpha_ohc_data$scenario <- "Hector - NMSE w/ unc (incl. OHC), Fit w/ S, Alpha" -hector_data <- rbind(default_data, nmse_ecs_data, nmse_bb_ecs_data, alpha_data, alpha_bb_data, alpha_vbb_data) +#hector_data <- rbind(default_data, nmse_data, nmse_bb_data, ecs_data, ecs_bb_data, alpha_data, alpha_bb_data) +hector_data <- rbind(default_data, alpha_data, alpha_ohc_data) hector_data$lower <- hector_data$value hector_data$upper <- hector_data$value diff --git a/scripts/initial_experiment.R b/scripts/initial_experiment.R index 41c9a93..5f35fd5 100644 --- a/scripts/initial_experiment.R +++ b/scripts/initial_experiment.R @@ -1,3 +1,4 @@ +# Script for experiment 1 # Script to re-implement the calibration of beta, q10, and diff from the V3 # manuscript # Author: Peter Scully @@ -22,7 +23,7 @@ TEMP_PATH <- INI_FILE <- system.file("input/hector_ssp245.ini", package = "hector") PARAMS <- c(BETA(), Q10_RH(), DIFFUSIVITY()) -OUTPUT <- file.path(RESULTS_DIR, "initial_experiment.txt") +OUTPUT <- file.path(RESULTS_DIR, "01_initial_experiment.txt") source(file.path(SCRIPTS_DIR, "major_functions.R")) diff --git a/scripts/initial_norm_exp.R b/scripts/initial_norm_exp.R index 9476097..7ea100d 100644 --- a/scripts/initial_norm_exp.R +++ b/scripts/initial_norm_exp.R @@ -1,3 +1,4 @@ +# Script for experiment 5 # Script to try normalizing the T and CO2 MSEs # manuscript # Author: Peter Scully @@ -22,7 +23,7 @@ TEMP_PATH <- INI_FILE <- system.file("input/hector_ssp245.ini", package = "hector") PARAMS <- c(BETA(), Q10_RH(), DIFFUSIVITY()) -OUTPUT <- file.path(RESULTS_DIR, "NMSE_big_box.txt") +OUTPUT <- file.path(RESULTS_DIR, "05_NMSE_big_box.txt") source(file.path(SCRIPTS_DIR, "major_functions.R")) diff --git a/scripts/major_functions.R b/scripts/major_functions.R index 779378c..629e60e 100644 --- a/scripts/major_functions.R +++ b/scripts/major_functions.R @@ -18,6 +18,9 @@ INPUT_TCR <- file.path(IDEAL_RUNS, "hector_1pctCO2.ini") OCEAN_AREA <- 5100656e8 * (1 - 0.29) # The total area of the ocean W_TO_ZJ <- 3.155693e-14 # Watts to ZJ +# Other Numeric Constants +CO2_OBS_SD <- 0.12 + # Constants for TCRE and ERF calculations EMISSIONS_VARS <- c(FFI_EMISSIONS(), LUC_EMISSIONS()) UPTAKE_VARS <- c(LUC_UPTAKE(), DACCS_UPTAKE()) @@ -64,8 +67,8 @@ get_co2_data <- function(file, scenario = "historical", include_unc = F) { co2_data$units <- " ppmv CO2" if (include_unc) { - co2_data$upper <- co2_data$value - co2_data$lower <- co2_data$value + co2_data$upper <- co2_data$value + CO2_OBS_SD + co2_data$lower <- co2_data$value - CO2_OBS_SD } return(co2_data) @@ -98,6 +101,49 @@ get_temp_data <- function(file, scenario = "historical", include_unc = F) { return(temp_data) } +# get_ohc_data - function to get ocean heat content data +# +# args: +# file - path to historical OHC data file +# scenario - name of scenario being run (default: "historical") +# include_unc - boolean indicating whether to include uncertainty data +# (default: F) +# +# returns: Hector-style data frame with OHC data +get_ohc_data <- function(file, scenario = "historical", include_unc = F) { + + # Reading in only OHC data + ohc_data <- read.table(file, + skip = 2, + sep = ",", + colClasses = c("numeric", "NULL", "NULL", "NULL", + "NULL", "NULL", "NULL", "numeric", + "numeric")) + + # Fixing table formatting + ohc_data <- na.omit(ohc_data) + colnames(ohc_data) <- c("year", "value", "unc") + + # Getting rid of non-integer years + ohc_data$year <- ohc_data$year - 0.5 + + # Adding in new columns to match Hector data frames + ohc_data$scenario <- scenario + ohc_data$variable <- "OHC" + ohc_data$units <- "ZJ" + + # Adding in confidence interval (if applicable) + if (include_unc) { + ohc_data$lower <- ohc_data$value - ohc_data$unc + ohc_data$upper <- ohc_data$value + ohc_data$unc + } + + # Getting rid of raw uncertainty column + ohc_data$unc <- NULL + + return(ohc_data) +} + ############################ ### Data Frame Functions ### @@ -257,6 +303,31 @@ run_hector <- function(ini_file, params, vals, yrs, vars, include_unc = F) { data <- rel_to_interval(data = data, var = GMST(), start = 1961, end = 1990) } + # Measuring ocean heat content (if applicable) + if (HEAT_FLUX() %in% vars) { + + ### Converting heat fluxes to cumulative heat content + + # Getting heat flux data + ohc_data <- filter(data, variable == HEAT_FLUX() & year %in% 1957:2014) + + # Converting heat flux to OHC change by year + ohc_data$value <- ohc_data$value * OCEAN_AREA * W_TO_ZJ + + # Converting OHC changes by year to total OHC (relative to 1957) + ohc_data$value <- cumsum(ohc_data$value) + ohc_data$variable <- "OHC" + + # Making OHC relative to 2005-2014 average + ohc_data <- rel_to_interval(data = ohc_data, + var = "OHC", + start = 2005, + end = 2014) + + # Adding ohc data to Hector data frame + data <- rbind(data, ohc_data) + } + # Adding in upper and lower bounds (if applicable) if (include_unc) { data$upper <- data$value diff --git a/scripts/more_params_unc_nmse.R b/scripts/more_params_unc_nmse.R index a292d41..93073de 100644 --- a/scripts/more_params_unc_nmse.R +++ b/scripts/more_params_unc_nmse.R @@ -1,3 +1,4 @@ +# Script for experiment 11 # Script to use normalized the T and CO2 MSEs while accounting for T uncertainty # Also includes ECS and alpha as params to optimize over # Author: Peter Scully @@ -22,7 +23,7 @@ TEMP_PATH <- INI_FILE <- system.file("input/hector_ssp245.ini", package = "hector") PARAMS <- c(BETA(), Q10_RH(), DIFFUSIVITY(), ECS(), AERO_SCALE()) -OUTPUT <- file.path(RESULTS_DIR, "alpha_ecs_unc_nmse_very_big_box.txt") +OUTPUT <- file.path(RESULTS_DIR, "11_alpha_ecs_unc_nmse_big_box.txt") source(file.path(SCRIPTS_DIR, "major_functions.R")) @@ -36,8 +37,8 @@ obs_data <- rbind(co2_data, temp_data) best_pars <- run_optim(obs_data = obs_data, ini_file = INI_FILE, params = PARAMS, - lower = c(0, 2.2 - 0.44 * 5, 2.3 - 0.1 * 5, 2 - 1, 0), - upper = c(0.5 + 0.232 * 5, 2.2 + 0.44 * 5, 2.3 + 0.1 * 5, 5 + 1, 3), + lower = c(0, 2.2 - 0.44 * 3, 2.3 - 0.1 * 3, 2, 0), + upper = c(0.5 + 0.232 * 3, 2.2 + 0.44 * 3, 2.3 + 0.1 * 3, 5, 3), yrs = 1750:2014, vars = c(GMST(), CONCENTRATIONS_CO2()), error_fn = mean_T_CO2_nmse_unc, diff --git a/scripts/other/ohc_comparison.R b/scripts/other/ohc_comparison.R index fa2e9de..a0f9fc2 100644 --- a/scripts/other/ohc_comparison.R +++ b/scripts/other/ohc_comparison.R @@ -44,6 +44,44 @@ source(file.path(SCRIPTS_DIR, "major_functions.R")) # 2005-2014. Note that new data frame will only contain OHC values, not # heat flux values calc_ohc <- function(params, vals, include_unc = F) { + + # Using an old version of run_hector for compatibility reasons + run_hector <- function(ini_file, params, vals, yrs, vars, include_unc = F) { + core <- newcore(ini_file) + + # Setting parameter values + if (!is.null(params)) { + for (i in 1:length(params)) { + setvar(core = core, + dates = NA, + var = params[i], + values = vals[i], + unit = getunits(params[i])) + } + } + reset(core) + + # Running core and fetching data + run(core) + data <- fetchvars(core, yrs, vars = vars) + shutdown(core) + + # Rescaling temperatures (if applicable) + if (GMST() %in% vars) { + data <- rel_to_interval(data = data, var = GMST(), start = 1961, end = 1990) + } + + # Adding in upper and lower bounds (if applicable) + if (include_unc) { + data$upper <- data$value + data$lower <- data$value + } + + return(data) + } + + + hect_data <- run_hector(INI_FILE, params = params, vals = vals, @@ -113,16 +151,16 @@ get_ohc_data <- function(file, scenario = "historical", include_unc = F) { obs_data <- get_ohc_data(OHC_FILE, include_unc = T) default_data <- calc_ohc(NULL, NULL, include_unc= T) -reg_box_data <- calc_ohc(PARAMS, c(0.732, 2.64, 2.2, 5, 1.15), include_unc= T) -big_box_data <- calc_ohc(PARAMS, c(1.196, 3.52, 2, 5, 0.948), include_unc= T) -bigger_box_data <- calc_ohc(PARAMS, c(1.66, 4.4, 1.8, 6, 0.83), include_unc= T) +reg_box_data <- calc_ohc(PARAMS, c(0.57, 1.76, 2.38, 2.95, 0.49), include_unc= T) +low_diff_data <- calc_ohc(PARAMS, c(0.57, 1.76, 1.1, 2.95, 0.49), include_unc= T) +ohc_optim_data <- calc_ohc(PARAMS, c(0.65, 1.76, 1.04, 2.33, 0.44), include_unc= T) default_data$scenario <- "Hector - Default" reg_box_data$scenario <- "Hector - All Params, Reg Box" -big_box_data$scenario <- "Hector - All Params, Big Box" -bigger_box_data$scenario <- "Hector - All Params, Bigger Box" +low_diff_data$scenario <- "Hector - All Params, Reg Box, Diff = 1.1" +ohc_optim_data$scenario <- "Hector - All Params, Reg Box, Matilda Diff, Optimize for OHC" -comb_data <- rbind(obs_data, default_data, reg_box_data, big_box_data, bigger_box_data) +comb_data <- rbind(obs_data, default_data, reg_box_data, low_diff_data, ohc_optim_data) ggplot(data = comb_data, aes(x = year, y = value, color = scenario)) + geom_ribbon(data = diff --git a/scripts/smoothed_temp_exp.R b/scripts/smoothed_temp_exp.R index 6aed236..3ad44f6 100644 --- a/scripts/smoothed_temp_exp.R +++ b/scripts/smoothed_temp_exp.R @@ -1,3 +1,4 @@ +# Scripts for experiments 2-4 # Script to run smoothed temperature runs # To keep running mean of a different number of observations, change ROLL_WIDTH # Author: Peter Scully @@ -27,7 +28,7 @@ PARAMS <- c(BETA(), Q10_RH(), DIFFUSIVITY()) ROLL_WIDTH <- 10 OUTPUT <- file.path(RESULTS_DIR, - paste("smooth_temps_", ROLL_WIDTH, "yr.txt", sep="")) + paste("04_smooth_temps_", ROLL_WIDTH, "yr.txt", sep="")) diff --git a/scripts/smoothed_temp_nmse.R b/scripts/smoothed_temp_nmse.R index 8f60f7c..a875c62 100644 --- a/scripts/smoothed_temp_nmse.R +++ b/scripts/smoothed_temp_nmse.R @@ -1,3 +1,4 @@ +# Script for experiments #6-#8 # Script to run smoothed temperature runs # To keep running mean of a different number of observations, change ROLL_WIDTH # Author: Peter Scully @@ -27,7 +28,7 @@ PARAMS <- c(BETA(), Q10_RH(), DIFFUSIVITY()) ROLL_WIDTH <- 10 OUTPUT <- file.path(RESULTS_DIR, - paste("smooth_temps_nmse_", ROLL_WIDTH, "yr_big_box.txt", sep="")) + paste("08_smooth_temps_nmse_", ROLL_WIDTH, "yr_big_box.txt", sep="")) @@ -52,7 +53,7 @@ obs_data <- rbind(co2_data, temp_data, smoothed_temp_data) best_pars <- run_optim(obs_data = obs_data, ini_file = INI_FILE, params = PARAMS, - lower = c(0, 2.2 - 0.44 * 3, 2.3 - 0.1 * 3), + lower = c(0, 2.2 - 0.44 * 3, 2.3 - 0.1 * 3), upper = c(0.5 + 0.232 * 3, 2.2 + 0.44 * 3, 2.3 + 0.1 * 3), yrs = 1750:2014, vars = c(GMST(), CONCENTRATIONS_CO2()), diff --git a/scripts/unc_nmse.R b/scripts/unc_nmse.R index 9f2db2e..7ae667d 100644 --- a/scripts/unc_nmse.R +++ b/scripts/unc_nmse.R @@ -1,3 +1,4 @@ +# Script for experiment #9 # Script to use normalized the T and CO2 MSEs while accounting for T uncertainty # manuscript # Author: Peter Scully @@ -22,7 +23,7 @@ TEMP_PATH <- INI_FILE <- system.file("input/hector_ssp245.ini", package = "hector") PARAMS <- c(BETA(), Q10_RH(), DIFFUSIVITY()) -OUTPUT <- file.path(RESULTS_DIR, "unc_nmse_big_box.txt") +OUTPUT <- file.path(RESULTS_DIR, "09_unc_nmse_big_box.txt") source(file.path(SCRIPTS_DIR, "major_functions.R")) @@ -36,8 +37,7 @@ obs_data <- rbind(co2_data, temp_data) best_pars <- run_optim(obs_data = obs_data, ini_file = INI_FILE, params = PARAMS, - par = c(0, 2.2 - 0.44*3, 2.6), - lower = c(0, 2.2 - 0.44 * 3, 2.3 - 0.1 * 3), + lower = c(0, 2.2 - 0.44 * 3, 2.3 - 0.1 * 3), upper = c(0.5 + 0.232 * 3, 2.2 + 0.44 * 3, 2.3 + 0.1 * 3), yrs = 1750:2014, vars = c(GMST(), CONCENTRATIONS_CO2()), diff --git a/tests/error_tests.R b/tests/error_tests.R index ce8b278..e075e5c 100644 --- a/tests/error_tests.R +++ b/tests/error_tests.R @@ -87,9 +87,21 @@ nmse_unc_tests <- function() { c(10, 7, 28)) == 9 / 625) } +# Trying to use mvsse +mvsse_tests <- function () { + # Just using mse tests and sd of 1 + assert_that(mvsse(0, 1, 0) == 0) + assert_that(mvsse(c(0, 1, 2), 1, c(1, 2, 1)) == 1) + assert_that(mvsse(c(0, 1, 2), c(1, 1, 1), c(0, 1, -1)) == 3) + assert_that(mvsse(c(-2, 1, 2), 1, c(0, 3, 0)) == 4) + + # Other tests + assert_that(mvsse(c(0, 1, 2, 3), c(2, 3, 4, 1), c(-1, 4, 4, 3)) == 0.375) +} # Calling all testing functions mse_tests() get_var_mse_tests() mse_unc_tests() -nmse_unc_tests() \ No newline at end of file +nmse_unc_tests() +mvsse_tests() \ No newline at end of file