Conformal bootstrap model

README.md

@@ -100,12 +100,27 @@ elexmodel 2017-11-07_VA_G --estimands=turnout --estimands my_estimand --estimand
 
 Some model types have specific model parameters that can be included.
 
-| Name      | Type    | Acceptable values           | model           |
-|-----------|---------|-----------------------------|-----------------|
-| lambda    | numeric | 0-inf                       | all             |
-| robust    | boolean | larger prediction intervals | `nonparametric` |
-| beta      | numeric | variance inflation          | `gaussian`      |
-| winsorize | boolean | winsorize std estimate      | `gaussian`      |
+| Name                              | Type    | Acceptable values                | model           |
+|-----------------------------------|---------|----------------------------------|-----------------|
+| lambda                            | numeric | regularization constant          | all             |
+| turnout_factor_lower              | numeric | drop units with < turnout factor | all             |
+| turnout_factor_upper              | numeric | drop units with < turnout factor | all             |
+| robust                            | boolean | larger prediction intervals      | `nonparametric` |
+| beta                              | numeric | variance inflation               | `gaussian`      |
+| winsorize                         | boolean | winsorize std estimate           | `gaussian`      |
+| B                                 | numeric | bootstrap samples                | `bootstrap`     |
+| T                                 | numeric | temperature for aggregate        | `bootstrap`     |
+| strata                            | list    | groups to stratify bootstrap by  | `bootstrap`     |
+| agg_model_hard_threshold          | bool    | hard threshold aggregate model   | `bootstrap`     |
+| y_LB                              | numeric | lower bound norm. margin dist    | `bootstrap`     |
+| y_UB                              | numeric | upper bound norm. margin dist    | `bootstrap`     |
+| z_LB                              | numeric | lower bound turnout fact dist    | `bootstrap`     |
+| z_UB                              | numeric | lower bound turnout fact dist    | `bootstrap`     |
+| y_unobserved_lower_bound          | numeric | lower bound for norm. margin     | `bootstrap`     |
+| y_unobserved_upper_bound          | numeric | upper bound for norm. margin     | `bootstrap`     |
+| percent_expected_vote_error_bound | numeric | error tolerance on expected vote | `bootstrap`     |
+| z_unobserved_lower_bound          | numeric | lower bound for turnout factor   | `bootstrap`     |
+| z_unobserved_upper_bound          | numeric | upper bound for turnout factor   | `bootstrap`     |
 
 This is the class and function that invokes the general function to generate estimates. You can install `elex-model` as a Python package and use this code snippet in other projects.
 

src/elexmodel/cli.py

@@ -35,12 +35,7 @@ def type_cast_value(self, ctx, value):
     "--pi_method",
     "pi_method",
     default="nonparametric",
-    type=click.Choice(
-        [
-            "gaussian",
-            "nonparametric",
-        ]
-    ),
+    type=click.Choice(["gaussian", "nonparametric", "bootstrap"]),
 )
 @click.option("--prediction_intervals", "prediction_intervals", default=[0.7, 0.9], multiple=True)
 @click.option("--percent_reporting_threshold", "percent_reporting_threshold", default=100)

src/elexmodel/client.py

@@ -10,6 +10,7 @@
 from elexmodel.handlers.data.ModelResults import ModelResultsHandler
 from elexmodel.handlers.data.PreprocessedData import PreprocessedDataHandler
 from elexmodel.logging import initialize_logging
+from elexmodel.models.BootstrapElectionModel import BootstrapElectionModel
 from elexmodel.models.ConformalElectionModel import ConformalElectionModel
 from elexmodel.models.GaussianElectionModel import GaussianElectionModel
 from elexmodel.models.NonparametricElectionModel import NonparametricElectionModel
@@ -39,6 +40,7 @@ def __init__(self):
         super().__init__()
         self.all_conformalization_data_unit_dict = defaultdict(dict)
         self.all_conformalization_data_agg_dict = defaultdict(dict)
+        self.model = None
 
     def _check_input_parameters(
         self,
@@ -88,8 +90,7 @@ def _check_input_parameters(
             ]
         if len(invalid_fixed_effects) > 0:
             raise ValueError(f"Fixed effect(s): {invalid_fixed_effects} not valid. Please check config")
-
-        if pi_method not in {"gaussian", "nonparametric"}:
+        if pi_method not in {"gaussian", "nonparametric", "bootstrap"}:
             raise ValueError(
                 f"Prediction interval method: {pi_method} is not valid. \
                     pi_method has to be either `gaussian` or `nonparametric`."
@@ -102,6 +103,14 @@ def _check_input_parameters(
                 not isinstance(model_parameters["lambda_"], (float, int)) or model_parameters["lambda_"] < 0
             ):
                 raise ValueError("lambda is not valid. It has to be numeric and greater than zero.")
+            if "turnout_factor_lower" in model_parameters and not isinstance(
+                model_parameters["turnout_factor_lower"], float
+            ):
+                raise ValueError("turnout_factor_lower is not valid. Has to be a float.")
+            if "turnout_factor_upper" in model_parameters and not isinstance(
+                model_parameters["turnout_factor_upper"], float
+            ):
+                raise ValueError("turnout_factor_upper is not valid. Has to be a float.")
             if pi_method == "gaussian":
                 if "beta" in model_parameters and not isinstance(model_parameters["beta"], (int, float)):
                     raise ValueError("beta is not valid. Has to be either an integer or a float.")
@@ -110,7 +119,45 @@ def _check_input_parameters(
             elif pi_method == "nonparametric":
                 if "robust" in model_parameters and not isinstance(model_parameters["robust"], bool):
                     raise ValueError("robust is not valid. Has to be a boolean.")
-
+            elif pi_method == "bootstrap":
+                if "B" in model_parameters and not isinstance(model_parameters["B"], int):
+                    raise ValueError("B is not valid. Has to be either an integer.")
+                if "T" in model_parameters and not isinstance(model_parameters["T"], (int, float)):
+                    raise ValueError("T is not valid. Has to be either an integer or a float.")
+                if "strata" in model_parameters and not isinstance(model_parameters["strata"], list):
+                    raise ValueError("strata is not valid. Has to be a list.")
+                if "agg_model_hard_threshold" in model_parameters and not isinstance(
+                    model_parameters["agg_model_hard_threshold"], bool
+                ):
+                    raise ValueError("agg_model_hard_threshold is not valid. Has to be a boolean.")
+                if "y_LB" in model_parameters and not isinstance(model_parameters["y_LB"], float):
+                    raise ValueError("y_LB is not valid. Has to be a float.")
+                if "y_UB" in model_parameters and not isinstance(model_parameters["y_UB"], float):
+                    raise ValueError("y_UB is not valid. Has to be a float.")
+                if "z_LB" in model_parameters and not isinstance(model_parameters["z_LB"], float):
+                    raise ValueError("z_LB is not valid. Has to be a float.")
+                if "z_UB" in model_parameters and not isinstance(model_parameters["z_UB"], float):
+                    raise ValueError("z_UB is not valid. Has to be a float.")
+                if "y_unobserved_upper_bound" in model_parameters and not isinstance(
+                    model_parameters["y_unobserved_upper_bound"], float
+                ):
+                    raise ValueError("y_unobserved_upper_bound is not valid. Has to be a float.")
+                if "y_unobserved_lower_bound" in model_parameters and not isinstance(
+                    model_parameters["y_unobserved_lower_bound"], float
+                ):
+                    raise ValueError("y_unobserved_lower_bound is not valid. Has to be a float.")
+                if "percent_expected_vote_error_bound" in model_parameters and not isinstance(
+                    model_parameters["percent_expected_vote_error_bound"], float
+                ):
+                    raise ValueError("z_UB is not valid. Has to be a float.")
+                if "z_unobserved_upper_bound" in model_parameters and not isinstance(
+                    model_parameters["z_unobserved_upper_bound"], float
+                ):
+                    raise ValueError("z_unobserved_upper_bound is not valid. Has to be a float.")
+                if "z_unobserved_lower_bound" in model_parameters and not isinstance(
+                    model_parameters["z_unobserved_lower_bound"], float
+                ):
+                    raise ValueError("z_unobserved_lower_bound is not valid. Has to be a float.")
         if handle_unreporting not in {"drop", "zero"}:
             raise ValueError("handle_unreporting must be either `drop` or `zero`")
 
@@ -122,6 +169,9 @@ def get_aggregate_list(self, office, aggregate):
         raw_aggregate_list = base_aggregate + [aggregate]
         return sorted(list(set(raw_aggregate_list)), key=lambda x: AGGREGATE_ORDER.index(x))
 
+    def get_national_summary_votes_estimates(self, nat_sum_data_dict=None, called_states={}, base_to_add=0, alpha=0.9):
+        return self.model.get_national_summary_estimates(nat_sum_data_dict, called_states, base_to_add, alpha=alpha)
+
     def get_estimates(
         self,
         current_data,  # list of lists
@@ -158,10 +208,15 @@ def get_estimates(
         save_conformalization = "conformalization" in save_output
         handle_unreporting = kwargs.get("handle_unreporting", "drop")
 
+        district_election = False
+        if office in {"H", "Y", "Z"}:
+            district_election = True
+
         model_settings = {
             "election_id": election_id,
             "office": office,
             "geographic_unit_type": geographic_unit_type,
+            "district_election": district_election,
             "features": features,
             "fixed_effects": fixed_effects,
             "save_conformalization": save_conformalization,
@@ -172,6 +227,7 @@ def get_estimates(
         config_handler = ConfigHandler(
             election_id, config=raw_config, s3_client=s3.S3JsonUtil(TARGET_BUCKET), save=save_config
         )
+
         self._check_input_parameters(
             config_handler,
             office,
@@ -184,6 +240,7 @@ def get_estimates(
             model_parameters,
             handle_unreporting,
         )
+
         states_with_election = config_handler.get_states(office)
         estimand_baselines = config_handler.get_estimand_baselines(office, estimands)
 
@@ -213,13 +270,18 @@ def get_estimates(
             handle_unreporting=handle_unreporting,
         )
 
+        turnout_factor_lower = model_parameters.get("turnout_factor_lower", 0.5)
+        turnout_factor_upper = model_parameters.get("turnout_factor_upper", 1.5)
+
         reporting_units = data.get_reporting_units(
-            percent_reporting_threshold, features_to_normalize=features, add_intercept=True
+            percent_reporting_threshold, turnout_factor_lower, turnout_factor_upper
         )
         nonreporting_units = data.get_nonreporting_units(
-            percent_reporting_threshold, features_to_normalize=features, add_intercept=True
+            percent_reporting_threshold, turnout_factor_lower, turnout_factor_upper
+        )
+        unexpected_units = data.get_unexpected_units(
+            percent_reporting_threshold, aggregates, turnout_factor_lower, turnout_factor_upper
         )
-        unexpected_units = data.get_unexpected_units(percent_reporting_threshold, aggregates)
 
         LOG.info(
             "Model parameters: \n prediction intervals: %s, percent reporting threshold: %s, \
@@ -232,13 +294,15 @@ def get_estimates(
         )
 
         if pi_method == "nonparametric":
-            model = NonparametricElectionModel(model_settings=model_settings)
+            self.model = NonparametricElectionModel(model_settings=model_settings)
         elif pi_method == "gaussian":
-            model = GaussianElectionModel(model_settings=model_settings)
+            self.model = GaussianElectionModel(model_settings=model_settings)
+        elif pi_method == "bootstrap":
+            self.model = BootstrapElectionModel(model_settings=model_settings)
 
         minimum_reporting_units_max = 0
         for alpha in prediction_intervals:
-            minimum_reporting_units = model.get_minimum_reporting_units(alpha)
+            minimum_reporting_units = self.model.get_minimum_reporting_units(alpha)
             if minimum_reporting_units > minimum_reporting_units_max:
                 minimum_reporting_units_max = minimum_reporting_units
 
@@ -270,24 +334,26 @@ def get_estimates(
         )
 
         for estimand in estimands:
-            unit_predictions = model.get_unit_predictions(reporting_units, nonreporting_units, estimand)
+            unit_predictions = self.model.get_unit_predictions(
+                reporting_units, nonreporting_units, estimand, unexpected_units=unexpected_units
+            )
             results_handler.add_unit_predictions(estimand, unit_predictions)
             # gets prediciton intervals for each alpha
             alpha_to_unit_prediction_intervals = {}
             for alpha in prediction_intervals:
-                alpha_to_unit_prediction_intervals[alpha] = model.get_unit_prediction_intervals(
+                alpha_to_unit_prediction_intervals[alpha] = self.model.get_unit_prediction_intervals(
                     results_handler.reporting_units, results_handler.nonreporting_units, alpha, estimand
                 )
-                if isinstance(model, ConformalElectionModel):
+                if isinstance(self.model, ConformalElectionModel):
                     self.all_conformalization_data_unit_dict[alpha][
                         estimand
-                    ] = model.get_all_conformalization_data_unit()
+                    ] = self.model.get_all_conformalization_data_unit()
 
             results_handler.add_unit_intervals(estimand, alpha_to_unit_prediction_intervals)
 
             for aggregate in results_handler.aggregates:
                 aggregate_list = self.get_aggregate_list(office, aggregate)
-                estimates_df = model.get_aggregate_predictions(
+                estimates_df = self.model.get_aggregate_predictions(
                     results_handler.reporting_units,
                     results_handler.nonreporting_units,
                     results_handler.unexpected_units,
@@ -296,7 +362,7 @@ def get_estimates(
                 )
                 alpha_to_agg_prediction_intervals = {}
                 for alpha in prediction_intervals:
-                    alpha_to_agg_prediction_intervals[alpha] = model.get_aggregate_prediction_intervals(
+                    alpha_to_agg_prediction_intervals[alpha] = self.model.get_aggregate_prediction_intervals(
                         results_handler.reporting_units,
                         results_handler.nonreporting_units,
                         results_handler.unexpected_units,
@@ -305,10 +371,10 @@ def get_estimates(
                         alpha_to_unit_prediction_intervals[alpha],
                         estimand,
                     )
-                    if isinstance(model, ConformalElectionModel):
+                    if isinstance(self.model, ConformalElectionModel):
                         self.all_conformalization_data_agg_dict[alpha][
                             estimand
-                        ] = model.get_all_conformalization_data_agg()
+                        ] = self.model.get_all_conformalization_data_agg()
 
                 # get all of the prediction intervals here
                 results_handler.add_agg_predictions(

src/elexmodel/handlers/config.py

@@ -69,11 +69,15 @@ def get_estimand_baselines(self, office, estimands):
         Return dict of baseline pointers for requested estimands
         """
         baseline_pointers = {estimand: self.get_baseline_pointer(office).get(estimand) for estimand in estimands}
+        if "margin" in estimands:
+            baseline_pointers["margin"] = "margin"
         return baseline_pointers
 
     def get_estimands(self, office):
         baseline_pointer = self.get_baseline_pointer(office)
         estimands = list(baseline_pointer.keys())
+        if self.election_id.endswith("G"):
+            estimands += ["margin"]  # would otherwise need to add margin to every single config
         return estimands
 
     def get_states(self, office):
@@ -92,7 +96,12 @@ def get_geographic_unit_types(self, office):
         return self._get_office_subconfig(office).get("geographic_unit_types")
 
     def get_features(self, office):
-        return self._get_office_subconfig(office).get("features", [])
+        features = self._get_office_subconfig(office).get("features", [])
+        if self.election_id.endswith("G"):
+            features += [
+                "baseline_normalized_margin"
+            ]  # would otherwise need to add baseline_margin to every single config
+        return features
 
     def get_aggregates(self, office):
         return self._get_office_subconfig(office).get("aggregates", [])