diff --git a/R/CalculationFunctions.R b/R/CalculationFunctions.R
index 30b450df..8139877d 100644
--- a/R/CalculationFunctions.R
+++ b/R/CalculationFunctions.R
@@ -93,10 +93,6 @@ prepareDemandVectorForImportResults <- function(model, demand = "Production", lo
     }
     # Calculate import demand vector y_m. 
     if(demand == "Production") {
-      # This option left in for validation purposes.
-      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") {
       y_m <- prepareImportConsumptionDemand(model, location = location)
diff --git a/R/DemandFunctions.R b/R/DemandFunctions.R
index 107baaa5..05cdd968 100644
--- a/R/DemandFunctions.R
+++ b/R/DemandFunctions.R
@@ -91,8 +91,7 @@ prepareDomesticProductionDemand <- function(model, location) {
 #' @return A named vector with demand
 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.")
+    y_m_p <- prepare2RDemand(model, location, domestic = FALSE, import = TRUE)
   } else {
     # Note that model$mu (i.e., ITA) is not included in import production demand
     # because it is included in Domestic Production Demand
diff --git a/R/StateiorFunctions.R b/R/StateiorFunctions.R
index 8e6721a4..970f2887 100644
--- a/R/StateiorFunctions.R
+++ b/R/StateiorFunctions.R
@@ -135,7 +135,7 @@ generate2RDirectRequirementsfromUseWithTrade <- function(model, domestic){
 }
 
 
-#' Prepares a production demand vector representing production for two region models
+#' Prepares a demand vector for two region models
 #' Demand for SoI = SoI2SoI + RoUS2SoI
 #' Demand for RoUS = SoI2RoUS + RoUS2RoUS
 #' @param model An EEIO model object with model specs and IO tables loaded
@@ -143,7 +143,7 @@ generate2RDirectRequirementsfromUseWithTrade <- function(model, domestic){
 #' @param domestic A logical parameter indicating whether to generate domestic demand vector.
 #' @param demand_type A str indicating whether demand is Production or Consumption
 #' @return A named vector with demand
-prepare2RDemand <- function(model, location, domestic, demand_type = "Production") {
+prepare2RDemand <- function(model, location, domestic, import=FALSE, demand_type = "Production") {
   # Get state abbreviations, e.g., "US-ME" and "RoUS"
   state_abb <- sub(".*/","",model$FinalDemandMeta$Code_Loc) ## Extract characters after /
   state_abb <- unique(state_abb)
@@ -154,9 +154,13 @@ prepare2RDemand <- function(model, location, domestic, demand_type = "Production
     # 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 if(import) {
+    # Subtracts domestic from total for each of 4 objects of the transactions lists
+    use_table <- Map(`-`, model$UseTransactionswithTrade, model$DomesticUseTransactionswithTrade)
   } else {
     use_table <- model$UseTransactionswithTrade
   }
+
   # Getting list of final demand columns used for the appropriate demand
   if(demand_type == "Production") {
     FD_columns <- unlist(sapply(list("HouseholdDemand", "InvestmentDemand",
@@ -176,22 +180,24 @@ prepare2RDemand <- function(model, location, domestic, demand_type = "Production
   # State with SoI 
   if(location == state_abb[1]) {
     # calculate final demand for SoI
-    if(domestic) {
-      SoI2SoI_y <- rowSums(use_table[["SoI2SoI"]][, c(FD_columns, ita_column, "ExportResidual")])
-    } else {
-      SoI2SoI_y <- rowSums(use_table[["SoI2SoI"]][, c(FD_columns, "ExportResidual")])
-    }
+    # if(domestic) {
+    #   SoI2SoI_y <- rowSums(use_table[["SoI2SoI"]][, c(FD_columns, ita_column, "ExportResidual")])
+    # } else {
+    #   SoI2SoI_y <- rowSums(use_table[["SoI2SoI"]][, c(FD_columns, "ExportResidual")])
+    # }
+    SoI2SoI_y <- rowSums(use_table[["SoI2SoI"]][, c(FD_columns, ita_column, "ExportResidual")])
     RoUS2SoI_y  <- rowSums(use_table[["RoUS2SoI"]][, c(FD_columns, ita_column)])
     y <- c(SoI2SoI_y, RoUS2SoI_y)
 
   # For RoUS
   } else if(location == state_abb[2]) {
     # calculate final demand for RoUS
-    if(domestic) {
-      RoUS2RoUS_y <- rowSums(use_table[["RoUS2RoUS"]][, c(FD_columns, ita_column, "ExportResidual")])
-    } else {
-      RoUS2RoUS_y <- rowSums(use_table[["RoUS2RoUS"]][, c(FD_columns, "ExportResidual")])
-    }
+    # if(domestic) {
+    #   RoUS2RoUS_y <- rowSums(use_table[["RoUS2RoUS"]][, c(FD_columns, ita_column, "ExportResidual")])
+    # } else {
+    #   RoUS2RoUS_y <- rowSums(use_table[["RoUS2RoUS"]][, c(FD_columns, "ExportResidual")])
+    # }
+    RoUS2RoUS_y <- rowSums(use_table[["RoUS2RoUS"]][, c(FD_columns, ita_column, "ExportResidual")])
     SoI2RoUS_y <- rowSums(use_table[["SoI2RoUS"]][, c(FD_columns, ita_column)])
     y <- c(SoI2RoUS_y, RoUS2RoUS_y)
   }
@@ -250,7 +256,7 @@ print2RValidationResults <- function(model) {
   twoRegionResults_ls$finalResultsProductionComplete <- finalResultsProductionComplete
   twoRegionResults_ls$finalResultsConsumptionComplete <- finalResultsConsumptionComplete
 
-  # return(twoRegionResults_ls)
+  return(twoRegionResults_ls)
 }