Merge branch 'develop' into two_region_disagg

.gitignore

@@ -7,3 +7,4 @@ work
 inst/doc/**/*.html
 renv/
 examples
+tests/*.xlsx

R/BuildModel.R

@@ -109,9 +109,11 @@ constructEEIOMatrices <- function(model, configpaths = NULL) {
   logging::loginfo("Calculating Rho matrix (price year ratio)...")
   model$Rho <- calculateModelIOYearbyYearPriceRatio(model)
 
-  # Calculate producer over purchaser price ratio.
-  logging::loginfo("Calculating Phi matrix (producer over purchaser price ratio)...")
-  model$Phi <- calculateProducerbyPurchaserPriceRatio(model)
+  if (model$specs$IODataSource!="stateior") { 
+    # Calculate producer over purchaser price ratio.
+    logging::loginfo("Calculating Phi matrix (producer over purchaser price ratio)...")
+    model$Phi <- calculateProducerbyPurchaserPriceRatio(model)
+  }
 
   # Calculate basic over producer price ratio.
   logging::loginfo("Calculating Tau matrix (basic over producer price ratio)...")
@@ -267,7 +269,7 @@ createCfromFactorsandBflows <- function(factors,B_flows) {
   C[is.na(C)] <- 0
 
   # Make sure CO2e flows are characterized (see issue #281)
-  f <- B_flows[!(B_flows %in% factors$Flow & grep("kg CO2e", B_flows))]  
+  f <- B_flows[!(B_flows %in% factors$Flow) & grepl("kg CO2e", B_flows)]
   C[, f] <- 1
   # Filter and resort model C flows and make it into a matrix
   C <- as.matrix(C[, B_flows])

R/CalculationFunctions.R

@@ -39,7 +39,8 @@ calculateEEIOModel <- function(model, perspective, demand = "Production", locati
 #' @param location, str optional location code for demand vector, required for two-region models
 #' @param use_domestic_requirements A logical value: if TRUE, use domestic demand and L_d matrix;
 #' if FALSE, use complete demand and L matrix.
-prepareDemandVectorForStandardResults <- function(model, demand = "Production", location = NULL, use_domestic_requirements = FALSE) {
+prepareDemandVectorForStandardResults <- function(model, demand = "Production",
+                                                  location = NULL, use_domestic_requirements = FALSE) {
   if (is.character(demand)) {
     #assume this is a model build-in demand 
     #try to load the model vector
@@ -90,11 +91,13 @@ prepareDemandVectorForImportResults <- function(model, demand = "Production", lo
       location <- "US"
     }
     # Calculate import demand vector y_m. 
-    if(demand == "Production"){
+    if(demand == "Production") {
       # This option left in for validation purposes.
-      logging::loginfo("Warning: Production demand vector not recommended for estimating results for models with external Import Factors. ")
+      logging::loginfo(paste0("Warning: Production demand vector not recommended ",
+                              "for estimating results for models with external ",
+                              "Import Factors."))
       y_m <- prepareImportProductionDemand(model, location = location)
-    } else if(demand == "Consumption"){
+    } else if(demand == "Consumption") {
       y_m <- prepareImportConsumptionDemand(model, location = location)
     }
   } else {
@@ -140,8 +143,7 @@ calculateResultsWithExternalFactors <- function(model, perspective = "FINAL", de
   rownames(hh_lcia) <- codes
 
   # Calculate Final perspective results
-  if(perspective == "FINAL"){
-
+  if(perspective == "FINAL") {
     y_d <- diag(as.vector(y_d))
     y_m <- diag(as.vector(y_m))
 
@@ -171,8 +173,7 @@ calculateResultsWithExternalFactors <- function(model, perspective = "FINAL", de
       result$LCIA_f <- rbind(result$LCIA_f, hh_lcia)
     }
 
-  } else{ # Calculate direct perspective results.
-
+  } else { # Calculate direct perspective results.
     # Direct perspective implemented using the following steps:
     # Imported_LCI = LCI_direct_domestic
     # Calculate Direct Perspective LCI (a matrix with total impacts in form of sector x flows)
@@ -186,7 +187,8 @@ calculateResultsWithExternalFactors <- function(model, perspective = "FINAL", de
       y_d <- diag(as.vector(y_d))
       y_m <- diag(as.vector(y_m))
 
-      totalLCI <- (model$B %*% model$L_d %*% y_d) + (model$Q_t %*% model$A_m %*% model$L_d %*% y_d + model$Q_t %*% y_m) # same equation as final perspective
+      totalLCI <- (model$B %*% model$L_d %*% y_d) + (model$Q_t %*% model$A_m %*% model$L_d %*% y_d + model$Q_t %*% y_m)
+      # ^^ same equation as final perspective
       totalLCI <- t(totalLCI)
       rownames(totalLCI) <- colnames(model$Q_t)
       # Taking the overall difference between domestic and total LCI, not by commodity, to be the totals for the imports.
@@ -201,8 +203,7 @@ calculateResultsWithExternalFactors <- function(model, perspective = "FINAL", de
       colnames(result$LCI_d) <- rownames(model$Q_t)
       rownames(result$LCI_d) <- colnames(model$Q_t)
     }
-
-
+
     # # Alternate approach for calculating Direct LCI for IF
     # # Attempting to duplicate calculateDirectPerspectiveLCI for IF models
     # econTermOne <- model$L_d %*% y_d
@@ -240,9 +241,8 @@ calculateResultsWithExternalFactors <- function(model, perspective = "FINAL", de
 }
 
 
-
-#' Calculate total emissions/resources (LCI) and total impacts (LCIA) for an EEIO model that does not have external import factors
-#' for a given perspective and demand vector.
+#' Calculate total emissions/resources (LCI) and total impacts (LCIA) for an EEIO model
+#' that does not have external import factors for a given perspective and demand vector.
 #' @param model A complete EEIO model: a list with USEEIO model components and attributes.
 #' @param perspective Perspective of the model, can be "DIRECT" or "FINAL". "DIRECT" perspective
 #' aligns results with the sectors in which they are produced, while "FINAL" perspective aligns
@@ -511,8 +511,12 @@ calculateSectorPurchasedbySectorSourcedImpact <- function (y, model, indicator)
 #' @return A list with M_margin and N_margin
 #' @export
 calculateMarginSectorImpacts <- function(model) {
+  if (model$specs$IODataSource == "stateior") {
+    stop("Margins not available for two-region models")
+  }
   # Calculation fractions of producer price for each margin
-  MarginCoefficients <- as.matrix(model$Margins[, c("Transportation", "Wholesale", "Retail")]/model$Margins[, c("ProducersValue")])
+  MarginCoefficients <- as.matrix(model$Margins[, c("Transportation", "Wholesale", "Retail")] / 
+                                    model$Margins[, c("ProducersValue")])
   rownames(MarginCoefficients) <- model$Margins$SectorCode
   MarginCoefficients[is.na(MarginCoefficients)] <- 0
 

R/DemandFunctions.R

@@ -92,9 +92,10 @@ prepareDomesticProductionDemand <- function(model, location) {
 prepareImportProductionDemand <- function(model, location) {
   if (model$specs$IODataSource == "stateior") {
     y_m_p <- prepare2RDemand(model, location, domestic = FALSE)
-     stop("Import production demand not yet implemented for 2R models.")
+    stop("Import production demand not yet implemented for 2R models.")
   } else {
-    # Note that model$mu (i.e., ITA) is not included in import production demand because it is included in Domestic Production Demand
+    # Note that model$mu (i.e., ITA) is not included in import production demand
+    # because it is included in Domestic Production Demand
     loc <- grepl(location, model$FinalDemandMeta$Code_Loc)
     export_code <- model$FinalDemandMeta[model$FinalDemandMeta$Group=="Export" & loc, "Code_Loc"]
     changeinventories_code <- model$FinalDemandMeta[model$FinalDemandMeta$Group=="ChangeInventories" & loc, "Code_Loc"]
@@ -130,12 +131,10 @@ prepareConsumptionDemand <- function(model, location) {
 #' @return a named vector with demand
 prepareImportConsumptionDemand <- function(model, location) {
   if (model$specs$IODataSource == "stateior") {
-    # y_c <- prepare2RDemand(model, location, domestic = FALSE, demand_type = "Consumption")
     ImportMatrix <- model$U - model$U_d
     ImportMatrix <- head(ImportMatrix, -6) # drop value add rows; TODO update this
     ImportFinalDemand <- ImportMatrix[, which(colnames(ImportMatrix) %in% model$FinalDemandMeta$Code_Loc)]
     y_c <- sumforConsumption(model, ImportFinalDemand, location)
-    # stop("Import consumption demand not yet implemented for 2R models.")
   } else {
     # Including InternationalTradeAdjustment in DomesticFinalDemand for import factors calculations
     ImportFinalDemand <- model$ImportMatrix[, which(colnames(model$ImportMatrix) %in% model$FinalDemandMeta$Code_Loc)]

R/DisaggregateFunctions.R

@@ -197,10 +197,14 @@ disaggregateSetup <- function (model, configpaths = NULL, setupType = "Disaggreg
     }
 
     # For Two-region model, develop two-region specs
-    if (model$specs$IODataSource=="stateior"){
+    if (model$specs$IODataSource=="stateior") {
       if (stringr::str_sub(spec$OriginalSectorCode, start=-3)=="/US") {
 
-        if(!is.null(spec$stateFile)){
+        # Create disaggregation specs from proxy data (e.g., employment by sector by state)
+        if(!is.null(spec$stateFile)){ 
+
+          stop("This section of code is meant to be used with 2R models with disaggregated utilities
+                           and is not yet fully implemented.")
 
           for(region in model$specs$ModelRegionAcronyms){
             # Define paramters for createDisaggFilesFromProxyData function call
@@ -240,7 +244,7 @@ disaggregateSetup <- function (model, configpaths = NULL, setupType = "Disaggreg
             specs[[regionDisaggCode]] <- d2
 
           }
-        } else{
+        } else {
           # Create disaggregation specs from national tables
           for(region in model$specs$ModelRegionAcronyms){
             d2 <- prepareTwoRegionDisaggregation(spec, region, model$specs$ModelRegionAcronyms)

R/LoadDemandVectors.R

@@ -7,15 +7,6 @@ loadDemandVectors <- function(model) {
   model$DemandVectors$vectors <- list()
   meta <- loadDefaultDemandVectorMeta(model)
 
-  # if(model$specs$ExternalImportFactors == TRUE){
-  #   ImportDemand <- list()
-  #   ImportDemand["Type"] <- "Production"
-  #   ImportDemand["System"] <- "Import"
-  #   
-  #   model$specs$DemandVectors$ImportDemandVectors <- ImportDemand
-  #   # Fails at the following line below: meta <- rbind(meta,data.frame(i, stringsAsFactors = FALSE) )
-  # }
-  # 
   specs <- model$specs$DemandVectors
   for (v in setdiff(names(specs), "DefaultDemand")) {
     # Populate metadata

R/LoadIOTables.R

@@ -41,8 +41,10 @@ loadIOData <- function(model, configpaths = NULL) {
     names(model$InternationalTradeAdjustment) <- model$Industries$Code_Loc
   }
 
-  # Add Margins table
-  model$Margins <- getMarginsTable(model)
+  if (model$specs$IODataSource != "stateior") {
+    # Add Margins table, currently only for one-region models (see Issue #290)
+    model$Margins <- getMarginsTable(model)
+  }
 
   # Add TaxLessSubsidies table
   model$TaxLessSubsidies <- generateTaxLessSubsidiesTable(model)

R/StateiorFunctions.R

@@ -139,7 +139,9 @@ prepare2RDemand <- function(model, location, domestic, demand_type = "Production
   state_abb <- unique(state_abb)
   iolevel <- model$specs$BaseIOLevel
 
-  if(domestic){#TODO: CHANGE domestic FROM BOOLEAN TO STRING WITH VALUES 'domestic', 'production', and 'import', so we can calculate the import matrix in the following if else if else block
+  if(domestic) {
+    # TODO: CHANGE domestic FROM BOOLEAN TO STRING WITH VALUES 'domestic', 'production', 
+    # and 'import', so we can calculate the import matrix in the following if else if else block
     use_table <- model$DomesticUseTransactionswithTrade    
   } else {
     use_table <- model$UseTransactionswithTrade
@@ -161,13 +163,6 @@ prepare2RDemand <- function(model, location, domestic, demand_type = "Production
       } else {
         SoI2SoI_y <- rowSums(use_table[["SoI2SoI"]][, c(FD_columns, "ExportResidual")])
       }
-
-      # if(!is.null(model$specs$ImportFactors)){
-      #   RoUS2SoI_y  <- rowSums(use_table[["RoUS2SoI"]][, c(FD_columns)]) # ITA column accounted for in domestic production demand vector for models with IF       
-      # }else{
-      #   RoUS2SoI_y  <- rowSums(use_table[["RoUS2SoI"]][, c(FD_columns, ita_column)])
-      # }
-      # 
       RoUS2SoI_y  <- rowSums(use_table[["RoUS2SoI"]][, c(FD_columns, ita_column)])
       y_p <- c(SoI2SoI_y, RoUS2SoI_y)
 
@@ -178,13 +173,6 @@ prepare2RDemand <- function(model, location, domestic, demand_type = "Production
       } else {
         RoUS2RoUS_y <- rowSums(use_table[["RoUS2RoUS"]][, c(FD_columns, "ExportResidual")])
       }
-
-      # if(!is.null(model$specs$ImportFactors)){
-      #   SoI2RoUS_y <- rowSums(use_table[["SoI2RoUS"]][, c(FD_columns)]) # ITA column accounted for in domestic production demand vector for models with IF       
-      # }else{
-      #   SoI2RoUS_y <- rowSums(use_table[["SoI2RoUS"]][, c(FD_columns, ita_column)])
-      # }
-      # 
       SoI2RoUS_y <- rowSums(use_table[["SoI2RoUS"]][, c(FD_columns, ita_column)])
       y_p <- c(SoI2RoUS_y, RoUS2RoUS_y)
     }
@@ -326,6 +314,12 @@ createDisaggFilesFromProxyData <- function(model, disagg, disaggYear, disaggStat
   #    in the Use table, the three Detail utility commodities (rows) will have that same split for across all columns (industries/final demand)
   # 2) The disagg parameter will contain a disagg$stateDF variable that includes the data for the relevant disaggState and disaggYear parameters.
 
+  if(!is.null(spec$stateFile)){ 
+
+    stop("This section of code is meant to be used with 2R models with disaggregated utilities
+                           and is not yet fully implemented.")
+  }
+
   #Get subset of ratios for current year
   stateDFYear <- subset(disagg$stateDF, Year == disaggYear & State == disaggState)
 
@@ -379,4 +373,4 @@ createDisaggFilesFromProxyData <- function(model, disagg, disaggYear, disaggStat
 
   return(model)
 
-}
+}

R/ValidateModel.R

@@ -474,7 +474,7 @@ validateImportFactorsApproach <- function(model, demand = "Consumption"){
 
   cat("\nTesting that LCIA results are equivalent between standard and coupled model approaches (i.e., LCIA = LCIA_dm) when\n")
   cat("assuming model$M = model$Q_t.\n")
-  all.equal(LCIA_dm, LCIA)
+  print(all.equal(LCIA_dm, LCIA))
 
 }
 

R/WriteModel.R

@@ -285,8 +285,11 @@ writeModelMetadata <- function(model, dirs) {
   years <- data.frame(ID=colnames(model$Rho), stringsAsFactors = FALSE)
   years$Index <- c(1:length(years$ID)-1)
   years <- years[, fields$years]
-  checkNamesandOrdering(years$ID, colnames(model$Phi),
-                        "years in years.csv and cols in Phi matrix")
+
+  if(!is.null(model$Phi)) {
+    checkNamesandOrdering(years$ID, colnames(model$Phi),
+                          "years in years.csv and cols in Phi matrix")
+  }
   checkNamesandOrdering(years$ID, colnames(model$Rho),
                         "years in years.csv and cols in Rho matrix")
   utils::write.csv(years, paste0(dirs$model, "/years.csv"), na = "",

tests/modelspecs/USEEIOv3.0-GHG.yml → tests/modelspecs/USEEIOv2.2-GHG.yml

@@ -1,4 +1,4 @@
-Model: "USEEIOv3.0-GHG" # 2017 Detail, Commodity
+Model: "USEEIOv2.2-GHG" # 2017 Detail, Commodity
 BaseIOSchema: 2017
 BaseIOLevel: "Detail"
 IOYear: 2017 # Year for IO data
@@ -43,17 +43,17 @@ Indicators:
     SimpleUnit: "Kilograms Carbon Dioxide (CO2)"
     SimpleName: "Greenhouse Gases"
     StaticSource: TRUE
-    StaticFile: "lciafmt/traci/TRACI_2.1_v1.0.0_5555779.parquet"
+    StaticFile: "lciafmt/ipcc/IPCC_v1.1.1_27ba917.parquet"
     FileLocation: "DataCommons"
     ScriptFunctionCall: "getImpactMethod" #function to call for script
     ScriptFunctionParameters: 
-      indicators: ["Global warming"]
+      indicators: ["AR5-100"]
     DataSources:
-      USEPA_TRACI_2.1:
-        Title: "TRACI 2.1"
-        Author: "USEPA"
-        DataYear: NA
-        URL: "https://www.epa.gov/chemical-research/tool-reduction-and-assessment-chemicals-and-other-environmental-impacts-traci"
+      IPCC_AR5:
+        Title: "IPCC Fifth Assessment Report: Direct Global Warming Potentials for 100 year time horizon"
+        Author: "IPCC"
+        DataYear: 2017
+        URL: ""
         Primary: TRUE
 
 DemandVectors:

tests/test_model_build.R

@@ -22,15 +22,15 @@ model <- useeior:::loadDemandVectors(model)
 model <- useeior:::constructEEIOMatrices(model)
 printValidationResults(model)
 
-## USEEIOv3.0-GHG Detail, commodity model (2017 Schema)
-m <- "USEEIOv3.0-GHG"
+## USEEIOv2.2-GHG Detail, commodity model (2017 Schema)
+m <- "USEEIOv2.2-GHG"
 cfg <- paste0("modelspecs/", m, ".yml")
 model <- buildModel(m, configpaths = file.path(cfg))
 printValidationResults(model)
 
-## USEEIOv3.0-s-GHG Summary, commodity model (2017 Schema)
+## USEEIOv2.2-s-GHG Summary, commodity model (2017 Schema)
 model <- useeior:::initializeModel(m, configpaths = file.path(cfg))
-model$specs$Model <- "USEEIOv3.0-s-GHG"
+model$specs$Model <- "USEEIOv2.2-s-GHG"
 model$specs$BaseIOLevel <- "Summary"
 model$crosswalk <- useeior:::getModelCrosswalk(model) # reassign for summary model
 model <- useeior:::loadIOData(model)
@@ -64,6 +64,7 @@ m <- "USEEIOv2.0-GHG"
 cfg <- paste0("modelspecs/", m, ".yml")
 model <- buildModel(m, configpaths = file.path(cfg))
 printValidationResults(model)
+writeModeltoXLSX(model, ".")
 
 ## USEEIOv2.0 Detail, industry model
 model <- useeior:::initializeModel(m, configpaths = file.path(cfg))
@@ -132,6 +133,7 @@ m <- "GAEEIOv1.0-s-GHG-19"
 cfg <- paste0("modelspecs/", m, ".yml")
 model <- buildModel(m, configpaths = file.path(cfg))
 printValidationResults(model)
+writeModeltoXLSX(model, ".")
 
 ## StateEEIOv1.0 Two-region Summary model with Import Factors
 cfg <- c(paste0("modelspecs/", m, ".yml"),