Make all damage variables the same sign

docs/src/science.md

@@ -216,12 +216,12 @@ The first part (rate) is always negative: As farmers have imperfect foresight an
 For the impact of the rate of climate change (i.e., the annual change of climate) on agriculture, the assumed model is:
 
 ```math
-A_{t,r}^{r} = \alpha_{r}\left( \frac{\Delta T_{t}}{0.04} \right)^{\beta} + \left( 1 - \frac{1}{\rho} \right)A_{t - 1,r}^{r} \tag{A.2}
+A_{t,r}^{r} = -\left(\alpha_{r}\left( \frac{\Delta T_{t}}{0.04} \right)^{\beta} + \left( \frac{1}{\rho} - 1 \right)A_{t - 1,r}^{r}\right) \tag{A.2}
 ```
 
 where
 
-  * ``A^{r}`` denotes damage in agricultural production as a fraction due the rate of climate change by time and region;
+  * ``A^{r}`` denotes damage in agricultural production as a fraction due the rate of climate change by time and region (a positive value means a decrease in production);
   * ``t`` denotes time;
   * ``r`` denotes region;
   * ``\Delta T`` denotes the change in the regional mean temperature (in degrees Celsius) between time ``t`` and ``t - 1``;
@@ -232,12 +232,12 @@ where
 The model for the impact due to the level of climate change since 1990 is:
 
 ```math
-A_{t,r}^{l} = \delta_{r}^{l}T_{t} + \delta_{r}^{q}T_{t}^{2} \tag{A.3}
+A_{t,r}^{l} = -\left(\delta_{r}^{l}T_{t} + \delta_{r}^{q}T_{t}^{2}\right) \tag{A.3}
 ```
 
 where
 
-  * ``A^{l}`` denotes the damage in agricultural production as a fraction due to the level of climate change by time and region;
+  * ``A^{l}`` denotes the damage in agricultural production as a fraction due to the level of climate change by time and region (a positive value means a decrease in production);
   * ``t`` denotes time;
   * ``r`` denotes region;
   * ``T`` denotes the change (in degree Celsius) in regional mean temperature relative to 1990;
@@ -246,12 +246,12 @@ where
 CO₂ fertilisation has a positive, but saturating effect on agriculture, specified by
 
 ```math
-A_{t,r}^{f} = \gamma_{r}\ln\frac{\text{CO}2_{t}}{275} \tag{A.4}
+A_{t,r}^{f} = -\left(\gamma_{r}\ln\frac{\text{CO}2_{t}}{275}\right) \tag{A.4}
 ```
 
 where
 
-  * ``A^{f}`` denotes damage in agricultural production as a fraction due to the CO2 fertilisation by time and region;
+  * ``A^{f}`` denotes damage in agricultural production as a fraction due to the CO2 fertilisation by time and region (a positive value means a decrease in production);
   * ``t`` denotes time;
   * ``r`` denotes region;
   * ``CO2`` denotes the atmospheric concentration of carbon dioxide (in parts per million by volume);
@@ -260,7 +260,7 @@ where
 
 The parameters in Table A are calibrated, following the procedure described in Tol (2002a), to the results of Kane *et al*. (1992), Reilly *et al*. (1994), Morita *et al*. (1994), Fischer *et al*. (1996), and Tsigas *et al*. (1996). These studies all use a global computable general equilibrium model, and report results with and without adaptation, and with and without CO₂ fertilisation. The regional results from these studies are assumed to hold for each country in the respective regions. They are averaged over the studies and the climate scenarios for each country, and aggregated to the *FUND* regions. The standard deviations in Table A follow from the spread between studies and scenarios. Equation (A.4) follows from the difference in results with and without CO2 fertilization. Equation (A.3) follows from the results with full adaptation. Equation (A.2) follows from the difference in results with and without adaptation.
 
-Equations (A.1)-(A.4) express the impact of climate change as a percentage of agricultural production. In order to express this as a percentage of income, we need to know the share of agricultural production in total income. This is assumed to fall with per capita income, that is,
+Equations (A.1)-(A.4) express the damages of climate change as a percentage of agricultural production. In order to express this as a percentage of income, we need to know the share of agricultural production in total income. This is assumed to fall with per capita income, that is,
 
 ```math
 \frac{\text{GA}P_{t,r}}{Y_{t,r}} = \frac{\text{GA}P_{1990,r}}{Y_{1990,r}}\left( \frac{y_{1990,r}}{y_{t,r}} \right)^{\epsilon} \tag{A.5}
@@ -284,12 +284,12 @@ The code for the agricultural impacts component can be found at [https://github.
 The model is:
 
 ```math
-F_{t,r} = \alpha_{r}\left( \frac{y_{t,r}}{y_{1990,r}} \right)^{\epsilon}\left( 0.5\left( \frac{T_{t}}{1.0} \right)^{\beta} + 0.5\gamma\ln\left( \frac{\text{CO}2_{t}}{275} \right) \right) \tag{F.1}
+F_{t,r} = -\alpha_{r}\left( \frac{y_{t,r}}{y_{1990,r}} \right)^{\epsilon}\left( 0.5\left( \frac{T_{t}}{1.0} \right)^{\beta} + 0.5\gamma\ln\left( \frac{\text{CO}2_{t}}{275} \right) \right) \tag{F.1}
 ```
 
 where
 
-  * ``F`` denotes the change in forestry consumer and producer surplus (as a share of total income);
+  * ``F`` denotes the damages in forestry consumer and producer surplus as a share of total income (a positive value means a decrease in surplus);
   * ``t`` denotes time;
   * ``r`` denotes region;
   * ``y`` denotes per capita income (in 1995 US dollar per person per year);
@@ -306,12 +306,12 @@ The parameter ``\alpha`` is estimated as the average of the estimates by Perez-G
 The impact of climate change on water resources follows:
 
 ```math
-W_{t,r} = \min\left\{ \alpha_{r}Y_{1990,r}\left( 1 - \tau \right)^{t - 2000}\left( \frac{y_{t,r}}{y_{1990,r}} \right)^{\beta}\left( \frac{P_{t,r}}{P_{1990,r}} \right)^{\eta}\left( \frac{T_{t}}{1.0} \right)^{\gamma},\frac{Y_{t,r}}{10} \right\} \tag{W.1}
+W_{t,r} = \max\left\{ -\alpha_{r}Y_{1990,r}\left( 1 - \tau \right)^{t - 2000}\left( \frac{y_{t,r}}{y_{1990,r}} \right)^{\beta}\left( \frac{P_{t,r}}{P_{1990,r}} \right)^{\eta}\left( \frac{T_{t}}{1.0} \right)^{\gamma},\frac{-Y_{t,r}}{10} \right\} \tag{W.1}
 ```
 
 where
 
-  * ``W`` denotes the change in water resources (in 1995 US dollar) at time ``t`` in region ``r``;
+  * ``W`` denotes the decrease in water resources (in 1995 US dollar) at time ``t`` in region ``r``;
   * ``t`` denotes time;
   * ``r`` denotes region;
   * ``y`` denotes per capita income (in 1995 US dollar) at time ``t`` in region ``r``;
@@ -330,12 +330,12 @@ These parameters are from calibrating *FUND* to the results of Downing *et al*.
 For space heating, the model is:
 
 ```math
-SH_{t,r} = \frac{\alpha_{r}Y_{1990,r}\frac{\operatorname{atan}T_{t}}{\operatorname{atan}{1.0}}\left( \frac{y_{t,r}}{y_{1990,r}} \right)^{\epsilon}\left( \frac{P_{t,r}}{P_{1990,r}} \right)}{\prod_{s = 1990}^{t}{\text{AEE}I_{s,r}}} \tag{E.1}
+SH_{t,r} = -\frac{\alpha_{r}Y_{1990,r}\frac{\operatorname{atan}T_{t}}{\operatorname{atan}{1.0}}\left( \frac{y_{t,r}}{y_{1990,r}} \right)^{\epsilon}\left( \frac{P_{t,r}}{P_{1990,r}} \right)}{\prod_{s = 1990}^{t}{\text{AEE}I_{s,r}}} \tag{E.1}
 ```
 
 where
 
-  * ``\text{SH}`` denotes the decrease in expenditure on space heating (in 1995 US dollar) at time ``t`` in region ``r``;
+  * ``\text{SH}`` denotes the increase in expenditure on space heating (in 1995 US dollar) at time ``t`` in region ``r``;
   * ``\text{t\ }``denotes time;
   * ``r`` denotes region;
   * ``Y`` denotes income (in 1995 US dollar) at time ``t`` in region ``r``;
@@ -351,7 +351,7 @@ These parameters are from calibrating *FUND* to the results of Downing *et al*.
 For space cooling, the model is:
 
 ```math
-SC_{t,r} = \frac{\alpha_{r}Y_{1990,r}\left( \frac{T_{t}}{1.0} \right)^{\beta}\left( \frac{y_{t,r}}{y_{1990,r}} \right)^{\epsilon}\left( \frac{P_{t,r}}{P_{1990,r}} \right)}{\prod_{s = 1990}^{t}{\text{AEE}I_{s,r}}} \tag{E.2}
+SC_{t,r} = -\frac{\alpha_{r}Y_{1990,r}\left( \frac{T_{t}}{1.0} \right)^{\beta}\left( \frac{y_{t,r}}{y_{1990,r}} \right)^{\epsilon}\left( \frac{P_{t,r}}{P_{1990,r}} \right)}{\prod_{s = 1990}^{t}{\text{AEE}I_{s,r}}} \tag{E.2}
 ```
 
 where

src/components/ImpactAggregationComponent.jl

@@ -38,11 +38,10 @@
         else
             for r in d.regions
                 v.eloss[t, r] = min(
-                    0.0 -
-                    p.water[t, r] -
-                    p.forests[t, r] -
-                    p.heating[t, r] -
-                    p.cooling[t, r] -
+                    p.water[t, r] +
+                    p.forests[t, r] +
+                    p.heating[t, r] +
+                    p.cooling[t, r] +
                     p.agcost[t, r] +
                     p.drycost[t, r] +
                     p.protcost[t, r] +

src/components/ImpactAgricultureComponent.jl

@@ -13,7 +13,7 @@
     agrate = Variable(index=[time,regions])
     aglevel = Variable(index=[time,regions])
     agco2 = Variable(index=[time,regions])
-    agcost = Variable(index=[time,regions])
+    agcost = Variable(index=[time,regions])  # Economic damages: positive value means damages, negative value means benefits 
 
     agrbm = Parameter(index=[regions])
     agtime = Parameter(index=[regions])
@@ -50,20 +50,20 @@
                 if isnan((dtemp / 0.04)^p.agnl)
                     v.agrate[t, r] = 0.0
                 else
-                    v.agrate[t, r] = p.agrbm[r] * (dtemp / 0.04)^p.agnl + (1.0 - 1.0 / p.agtime[r]) * v.agrate[t - 1, r]
+                    v.agrate[t, r] = -1 * (p.agrbm[r] * (dtemp / 0.04)^p.agnl + (1.0 / p.agtime[r] - 1.0) * v.agrate[t - 1, r])
                 end
             end
 
             for r in d.regions
-                v.aglevel[t, r] = p.aglparl[r] * p.temp[t, r] + p.aglparq[r] * p.temp[t, r]^2.0
+                v.aglevel[t, r] = -1 * (p.aglparl[r] * p.temp[t, r] + p.aglparq[r] * p.temp[t, r]^2.0)
             end
 
             for r in d.regions
-                v.agco2[t, r] = p.agcbm[r] / log(2.0) * log(p.acco2[t - 1] / p.co2pre)
+                v.agco2[t, r] = -1 * p.agcbm[r] / log(2.0) * log(p.acco2[t - 1] / p.co2pre)
             end
 
             for r in d.regions
-                v.agcost[t, r] = min(1.0, v.agrate[t, r] + v.aglevel[t, r] + v.agco2[t, r]) * v.agrish[t, r] * p.income[t, r]
+                v.agcost[t, r] = max(-1.0, v.agrate[t, r] + v.aglevel[t, r] + v.agco2[t, r]) * v.agrish[t, r] * p.income[t, r]
             end
         end
     end

src/components/ImpactCoolingComponent.jl

@@ -1,7 +1,7 @@
 @defcomp impactcooling begin
     regions = Index()
 
-    cooling = Variable(index=[time,regions])
+    cooling = Variable(index=[time,regions])    # Economic damages: positive value means an increase in expenditures, negative value means a decrease
 
     cebm = Parameter(index=[regions])
     gdp90 = Parameter(index=[regions])
@@ -24,7 +24,7 @@
                 ypc = p.income[t, r] / p.population[t, r] * 1000.0
                 ypc90 = p.gdp90[r] / p.pop90[r] * 1000.0
 
-                v.cooling[t, r] = p.cebm[r] * p.cumaeei[t, r] * p.gdp90[r] * (p.temp[t, r] / 1.0)^p.cenl * (ypc / ypc90)^p.ceel * p.population[t, r] / p.pop90[r]
+                v.cooling[t, r] = -1 * p.cebm[r] * p.cumaeei[t, r] * p.gdp90[r] * (p.temp[t, r] / 1.0)^p.cenl * (ypc / ypc90)^p.ceel * p.population[t, r] / p.pop90[r]
             end
         end
     end

src/components/ImpactForestsComponent.jl

@@ -1,7 +1,7 @@
 @defcomp impactforests begin
     regions = Index()
 
-    forests = Variable(index=[time,regions])
+    forests = Variable(index=[time,regions])    # Economic damages: positive value means damages, negative value means benefits
 
     forbm = Parameter(index=[regions])
     gdp90 = Parameter(index=[regions])
@@ -26,11 +26,10 @@
                 ypc90 = p.gdp90[r] / p.pop90[r] * 1000.0
 
                 # TODO -oDavid Anthoff: RT uses -lP.forel for ppp case
-                v.forests[t, r] = p.forbm[r] * p.income[t, r] * (ypc / ypc90)^p.forel * (0.5 * p.temp[t, r]^p.fornl + 0.5 * log(p.acco2[t - 1] / p.co2pre) * p.forco2)
-
-                if v.forests[t, r] > 0.1 * p.income[t, r]
-                    v.forests[t, r] = 0.1 * p.income[t, r]
-                end
+                v.forests[t, r] = max(
+                    -1 * p.forbm[r] * p.income[t, r] * (ypc / ypc90)^p.forel * (0.5 * p.temp[t, r]^p.fornl + 0.5 * log(p.acco2[t - 1] / p.co2pre) * p.forco2),
+                    -0.1 * p.income[t, r]
+                )
             end
         end
     end

src/components/ImpactHeatingComponent.jl

@@ -1,7 +1,7 @@
 @defcomp impactheating begin
     regions = Index()
 
-    heating = Variable(index=[time,regions])
+    heating = Variable(index=[time,regions])    # Economic damages: positive value means an increase in expenditures, negative value means a decrease
 
     hebm = Parameter(index=[regions])
     gdp90 = Parameter(index=[regions])
@@ -24,7 +24,7 @@
                 ypc = p.income[t, r] / p.population[t, r] * 1000.0
                 ypc90 = p.gdp90[r] / p.pop90[r] * 1000.0
 
-                v.heating[t, r] = p.hebm[r] * p.cumaeei[t, r] * p.gdp90[r] * atan(p.temp[t, r]) / atan(1.0) * (ypc / ypc90)^p.heel * p.population[t, r] / p.pop90[r]
+                v.heating[t, r] = -1 * p.hebm[r] * p.cumaeei[t, r] * p.gdp90[r] * atan(p.temp[t, r]) / atan(1.0) * (ypc / ypc90)^p.heel * p.population[t, r] / p.pop90[r]
             end
         end
     end

src/components/ImpactWaterResourcesComponent.jl

@@ -4,7 +4,7 @@
     watech = Variable(index=[time])
     watechrate = Parameter(default = 0.005)
 
-    water = Variable(index=[time,regions])
+    water = Variable(index=[time,regions])  # Economic damages: positive value means damages, negative value means benefits
     wrbm = Parameter(index=[regions])
     wrel = Parameter(default = 0.85)
     wrnl = Parameter(default = 1)
@@ -30,13 +30,10 @@
             ypc = p.income[t, r] / p.population[t, r] * 1000.0
             ypc90 = p.gdp90[r] / p.pop90[r] * 1000.0
 
-            water = p.wrbm[r] * p.gdp90[r] * v.watech[t] * (ypc / ypc90)^p.wrel * (p.population[t, r] / p.pop90[r])^p.wrpl * p.temp[t, r]^p.wrnl
-
-            if water > 0.1 * p.income[t, r]
-                v.water[t, r] = 0.1 * p.income[t, r]
-            else
-                v.water[t, r] = water
-            end
+            v.water[t, r] = max(
+                -1 * p.wrbm[r] * p.gdp90[r] * v.watech[t] * (ypc / ypc90)^p.wrel * (p.population[t, r] / p.pop90[r])^p.wrpl * p.temp[t, r]^p.wrnl,
+                -0.1 * p.income[t, r]
+            )
         end
     end
 end

test/runtests.jl

@@ -13,16 +13,10 @@ using Query
 
 @testset "fund-model" begin
 
-#default model exported by fund module
-default_nsteps = 1050
-m = MimiFUND.get_model()
-run(m)
-@test Mimi.time_labels(m) == collect(1950:1:1950+default_nsteps)
-
 #default model created by MimiFUND.get_model()
 m1 = MimiFUND.get_model()
 run(m1)
-@test Mimi.time_labels(m1) == collect(1950:1:1950+default_nsteps)
+@test Mimi.time_labels(m1) == collect(1950:3000)
 
 # use optional `nsteps` arg for MimiFUND.get_model()
 @test_throws ErrorException m2 = MimiFUND.get_model(nsteps = 2000) # should error because it's longer than the default without providing different parameters
@@ -40,6 +34,10 @@ end #fund-model testset
 #------------------------------------------------------------------------------
 @testset "test-integration" begin
 
+# The signs of the following variables have been switched since the validation data were saved;
+#   need to multiply by -1 in their respective tests
+_updated_vars = [:agcost, :cooling, :heating, :forests, :water]
+
 m = MimiFUND.get_model()
 run(m)
 
@@ -53,7 +51,7 @@ for c in Mimi.compdefs(m), v in Mimi.variable_names(m, nameof(c))
     # load data for comparison
     orig_name = c.comp_id.comp_name
     filename = joinpath(@__DIR__, "../contrib/validation_data_v040/$orig_name-$v.csv")
-    results = m[nameof(c), v]
+    results = v in _updated_vars ? -1 * m[nameof(c), v] : m[nameof(c), v]
 
     df = load(filename) |> DataFrame
     if typeof(results) <: Number