Created
June 8, 2023 06:32
-
-
Save lmassaron/a01430584d001b8513661d9b51fe9527 to your computer and use it in GitHub Desktop.
Poisson regression
This file contains hidden or bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
| import math | |
| class RunningStats: | |
| def __init__(self): | |
| self.count = 0 | |
| self.mean = 0.0 | |
| self.M2 = 0.0 | |
| def update(self, value): | |
| self.count += 1 | |
| delta = value - self.mean | |
| self.mean += delta / self.count | |
| delta2 = value - self.mean | |
| self.M2 += delta * delta2 | |
| def get_mean(self): | |
| return self.mean | |
| def get_variance(self): | |
| if self.count < 2: | |
| return float('nan') | |
| else: | |
| return self.M2 / (self.count - 1) | |
| def get_standard_deviation(self): | |
| return math.sqrt(self.get_variance()) | |
| from sklearn.linear_model import PoissonRegressor | |
| from sklearn.utils import resample | |
| # Fit the Poisson regression model | |
| X = ... # Your design matrix | |
| y = ... # Your response variable | |
| poisson_model = PoissonRegressor() | |
| poisson_model.fit(X, y) | |
| # Generate predictions | |
| X_new = ... # New data for prediction | |
| n_bootstrap = 1000 # Number of bootstrap samples | |
| bootstrap_predictions = [] | |
| for _ in range(n_bootstrap): | |
| # Resample the data | |
| X_boot, y_boot = resample(X, y) | |
| # Fit the model to the resampled data | |
| boot_model = PoissonRegressor() | |
| boot_model.fit(X_boot, y_boot) | |
| # Make predictions on the new data | |
| boot_predictions = boot_model.predict(X_new) | |
| bootstrap_predictions.append(boot_predictions) | |
| # Calculate confidence intervals | |
| confidence_level = 0.95 # Specify the desired confidence level | |
| alpha = 1 - confidence_level | |
| lower_percentile = (alpha / 2) * 100 | |
| upper_percentile = (1 - alpha / 2) * 100 | |
| lower_bound = np.percentile(bootstrap_predictions, lower_percentile, axis=0) | |
| upper_bound = np.percentile(bootstrap_predictions, upper_percentile, axis=0) | |
| # Print the confidence intervals | |
| for i, (lower, upper) in enumerate(zip(lower_bound, upper_bound)): | |
| print(f"Prediction {i+1}: [{lower}, {upper}]") |
Sign up for free
to join this conversation on GitHub.
Already have an account?
Sign in to comment