Skip to content

Commit 3acc970

Browse files
cetagostinicetagostini
committed
Add notebook copy file
1 parent af4864a commit 3acc970

File tree

1 file changed

+52
-10
lines changed

1 file changed

+52
-10
lines changed

‎docs/source/notebooks/its_pymc copy.ipynb‎

Lines changed: 52 additions & 10 deletions
Original file line numberDiff line numberDiff line change
@@ -169,22 +169,64 @@
169169
},
170170
{
171171
"cell_type": "code",
172-
"execution_count": 6,
172+
"execution_count": 4,
173173
"metadata": {},
174174
"outputs": [
175175
{
176-
"ename": "ShapeError",
177-
"evalue": "Length of `dims` must match the dimensions of the dataset. (actual 1 != expected 2)",
176+
"name": "stderr",
177+
"output_type": "stream",
178+
"text": [
179+
"Initializing NUTS using jitter+adapt_diag...\n",
180+
"Multiprocess sampling (4 chains in 4 jobs)\n",
181+
"NUTS: [beta, sigma]\n"
182+
]
183+
},
184+
{
185+
"data": {
186+
"application/vnd.jupyter.widget-view+json": {
187+
"model_id": "c7484a104e064b2faf0bacbd39297b2c",
188+
"version_major": 2,
189+
"version_minor": 0
190+
},
191+
"text/plain": [
192+
"Output()"
193+
]
194+
},
195+
"metadata": {},
196+
"output_type": "display_data"
197+
},
198+
{
199+
"data": {
200+
"text/html": [
201+
"<pre style=\"white-space:pre;overflow-x:auto;line-height:normal;font-family:Menlo,'DejaVu Sans Mono',consolas,'Courier New',monospace\"></pre>\n"
202+
],
203+
"text/plain": []
204+
},
205+
"metadata": {},
206+
"output_type": "display_data"
207+
},
208+
{
209+
"name": "stderr",
210+
"output_type": "stream",
211+
"text": [
212+
"Sampling 4 chains for 1_000 tune and 1_000 draw iterations (4_000 + 4_000 draws total) took 2 seconds.\n",
213+
"The rhat statistic is larger than 1.01 for some parameters. This indicates problems during sampling. See https://arxiv.org/abs/1903.08008 for details\n",
214+
"Sampling: [beta, sigma, y_hat]\n",
215+
"Sampling: [y_hat]\n",
216+
"Sampling: [y_hat]\n"
217+
]
218+
},
219+
{
220+
"ename": "AttributeError",
221+
"evalue": "'numpy.ndarray' object has no attribute 'sel'",
178222
"output_type": "error",
179223
"traceback": [
180224
"\u001b[31m---------------------------------------------------------------------------\u001b[39m",
181-
"\u001b[31mShapeError\u001b[39m Traceback (most recent call last)",
182-
"\u001b[36mCell\u001b[39m\u001b[36m \u001b[39m\u001b[32mIn[6]\u001b[39m\u001b[32m, line 1\u001b[39m\n\u001b[32m----> \u001b[39m\u001b[32m1\u001b[39m result = \u001b[43mcp\u001b[49m\u001b[43m.\u001b[49m\u001b[43mStructuralTimeSeries\u001b[49m\u001b[43m(\u001b[49m\n\u001b[32m 2\u001b[39m \u001b[43m \u001b[49m\u001b[43mdf\u001b[49m\u001b[43m,\u001b[49m\n\u001b[32m 3\u001b[39m \u001b[43m \u001b[49m\u001b[43mtreatment_time\u001b[49m\u001b[43m,\u001b[49m\n\u001b[32m 4\u001b[39m \u001b[43m \u001b[49m\u001b[43mformula\u001b[49m\u001b[43m=\u001b[49m\u001b[33;43m\"\u001b[39;49m\u001b[33;43my ~ 1 + t + C(month)\u001b[39;49m\u001b[33;43m\"\u001b[39;49m\u001b[43m,\u001b[49m\n\u001b[32m 5\u001b[39m \u001b[43m \u001b[49m\u001b[43mmodel\u001b[49m\u001b[43m=\u001b[49m\u001b[43mcp\u001b[49m\u001b[43m.\u001b[49m\u001b[43mpymc_models\u001b[49m\u001b[43m.\u001b[49m\u001b[43mLinearRegression\u001b[49m\u001b[43m(\u001b[49m\u001b[43msample_kwargs\u001b[49m\u001b[43m=\u001b[49m\u001b[43m{\u001b[49m\u001b[33;43m\"\u001b[39;49m\u001b[33;43mrandom_seed\u001b[39;49m\u001b[33;43m\"\u001b[39;49m\u001b[43m:\u001b[49m\u001b[43m \u001b[49m\u001b[43mseed\u001b[49m\u001b[43m}\u001b[49m\u001b[43m)\u001b[49m\u001b[43m,\u001b[49m\n\u001b[32m 6\u001b[39m \u001b[43m)\u001b[49m\n",
183-
"\u001b[36mFile \u001b[39m\u001b[32m~/Documents/GitHub/CausalPy/causalpy/experiments/structural_time_series.py:129\u001b[39m, in \u001b[36mStructuralTimeSeries.__init__\u001b[39m\u001b[34m(self, data, treatment_time, formula, model, **kwargs)\u001b[39m\n\u001b[32m 124\u001b[39m \u001b[38;5;66;03m# Determine X to pass: if pre_X has 0 columns (e.g. y ~ 1), pass None or empty array.\u001b[39;00m\n\u001b[32m 125\u001b[39m \u001b[38;5;66;03m# This depends on how BSTS's _prepare_time_and_exog_features handles it.\u001b[39;00m\n\u001b[32m 126\u001b[39m \u001b[38;5;66;03m# Assuming it handles None or 0-column array gracefully if self.labels is empty.\u001b[39;00m\n\u001b[32m 127\u001b[39m X_fit = \u001b[38;5;28mself\u001b[39m.pre_X \u001b[38;5;28;01mif\u001b[39;00m \u001b[38;5;28mself\u001b[39m.pre_X.shape[\u001b[32m1\u001b[39m] > \u001b[32m0\u001b[39m \u001b[38;5;28;01melse\u001b[39;00m \u001b[38;5;28;01mNone\u001b[39;00m\n\u001b[32m--> \u001b[39m\u001b[32m129\u001b[39m \u001b[38;5;28;43mself\u001b[39;49m\u001b[43m.\u001b[49m\u001b[43mmodel\u001b[49m\u001b[43m.\u001b[49m\u001b[43mfit\u001b[49m\u001b[43m(\u001b[49m\u001b[43mX\u001b[49m\u001b[43m=\u001b[49m\u001b[43mX_fit\u001b[49m\u001b[43m,\u001b[49m\u001b[43m \u001b[49m\u001b[43my\u001b[49m\u001b[43m=\u001b[49m\u001b[38;5;28;43mself\u001b[39;49m\u001b[43m.\u001b[49m\u001b[43mpre_y\u001b[49m\u001b[43m,\u001b[49m\u001b[43m \u001b[49m\u001b[43mcoords\u001b[49m\u001b[43m=\u001b[49m\u001b[43mpre_coords\u001b[49m\u001b[43m)\u001b[49m\n\u001b[32m 131\u001b[39m \u001b[38;5;28;01melif\u001b[39;00m \u001b[38;5;28misinstance\u001b[39m(\u001b[38;5;28mself\u001b[39m.model, RegressorMixin):\n\u001b[32m 132\u001b[39m \u001b[38;5;28mself\u001b[39m.model.fit(X=\u001b[38;5;28mself\u001b[39m.pre_X, y=\u001b[38;5;28mself\u001b[39m.pre_y)\n",
184-
"\u001b[36mFile \u001b[39m\u001b[32m~/Documents/GitHub/CausalPy/causalpy/pymc_models.py:114\u001b[39m, in \u001b[36mPyMCModel.fit\u001b[39m\u001b[34m(self, X, y, coords)\u001b[39m\n\u001b[32m 110\u001b[39m \u001b[38;5;66;03m# Ensure random_seed is used in sample_prior_predictive() and\u001b[39;00m\n\u001b[32m 111\u001b[39m \u001b[38;5;66;03m# sample_posterior_predictive() if provided in sample_kwargs.\u001b[39;00m\n\u001b[32m 112\u001b[39m random_seed = \u001b[38;5;28mself\u001b[39m.sample_kwargs.get(\u001b[33m\"\u001b[39m\u001b[33mrandom_seed\u001b[39m\u001b[33m\"\u001b[39m, \u001b[38;5;28;01mNone\u001b[39;00m)\n\u001b[32m--> \u001b[39m\u001b[32m114\u001b[39m \u001b[38;5;28;43mself\u001b[39;49m\u001b[43m.\u001b[49m\u001b[43mbuild_model\u001b[49m\u001b[43m(\u001b[49m\u001b[43mX\u001b[49m\u001b[43m,\u001b[49m\u001b[43m \u001b[49m\u001b[43my\u001b[49m\u001b[43m,\u001b[49m\u001b[43m \u001b[49m\u001b[43mcoords\u001b[49m\u001b[43m)\u001b[49m\n\u001b[32m 115\u001b[39m \u001b[38;5;28;01mwith\u001b[39;00m \u001b[38;5;28mself\u001b[39m:\n\u001b[32m 116\u001b[39m \u001b[38;5;28mself\u001b[39m.idata = pm.sample(**\u001b[38;5;28mself\u001b[39m.sample_kwargs)\n",
185-
"\u001b[36mFile \u001b[39m\u001b[32m~/Documents/GitHub/CausalPy/causalpy/pymc_models.py:257\u001b[39m, in \u001b[36mLinearRegression.build_model\u001b[39m\u001b[34m(self, X, y, coords)\u001b[39m\n\u001b[32m 255\u001b[39m \u001b[38;5;28mself\u001b[39m.add_coords(coords)\n\u001b[32m 256\u001b[39m X = pm.Data(\u001b[33m\"\u001b[39m\u001b[33mX\u001b[39m\u001b[33m\"\u001b[39m, X, dims=[\u001b[33m\"\u001b[39m\u001b[33mobs_ind\u001b[39m\u001b[33m\"\u001b[39m, \u001b[33m\"\u001b[39m\u001b[33mcoeffs\u001b[39m\u001b[33m\"\u001b[39m])\n\u001b[32m--> \u001b[39m\u001b[32m257\u001b[39m y = \u001b[43mpm\u001b[49m\u001b[43m.\u001b[49m\u001b[43mData\u001b[49m\u001b[43m(\u001b[49m\u001b[33;43m\"\u001b[39;49m\u001b[33;43my\u001b[39;49m\u001b[33;43m\"\u001b[39;49m\u001b[43m,\u001b[49m\u001b[43m \u001b[49m\u001b[43my\u001b[49m\u001b[43m,\u001b[49m\u001b[43m \u001b[49m\u001b[43mdims\u001b[49m\u001b[43m=\u001b[49m\u001b[33;43m\"\u001b[39;49m\u001b[33;43mobs_ind\u001b[39;49m\u001b[33;43m\"\u001b[39;49m\u001b[43m)\u001b[49m\n\u001b[32m 258\u001b[39m beta = pm.Normal(\u001b[33m\"\u001b[39m\u001b[33mbeta\u001b[39m\u001b[33m\"\u001b[39m, \u001b[32m0\u001b[39m, \u001b[32m50\u001b[39m, dims=\u001b[33m\"\u001b[39m\u001b[33mcoeffs\u001b[39m\u001b[33m\"\u001b[39m)\n\u001b[32m 259\u001b[39m sigma = pm.HalfNormal(\u001b[33m\"\u001b[39m\u001b[33msigma\u001b[39m\u001b[33m\"\u001b[39m, \u001b[32m1\u001b[39m)\n",
186-
"\u001b[36mFile \u001b[39m\u001b[32m/opt/anaconda3/envs/CausalPy/lib/python3.13/site-packages/pymc/data.py:393\u001b[39m, in \u001b[36mData\u001b[39m\u001b[34m(name, value, dims, coords, infer_dims_and_coords, mutable, **kwargs)\u001b[39m\n\u001b[32m 391\u001b[39m dims = (dims,)\n\u001b[32m 392\u001b[39m \u001b[38;5;28;01mif\u001b[39;00m \u001b[38;5;129;01mnot\u001b[39;00m (dims \u001b[38;5;129;01mis\u001b[39;00m \u001b[38;5;28;01mNone\u001b[39;00m \u001b[38;5;129;01mor\u001b[39;00m \u001b[38;5;28mlen\u001b[39m(dims) == x.ndim):\n\u001b[32m--> \u001b[39m\u001b[32m393\u001b[39m \u001b[38;5;28;01mraise\u001b[39;00m pm.exceptions.ShapeError(\n\u001b[32m 394\u001b[39m \u001b[33m\"\u001b[39m\u001b[33mLength of `dims` must match the dimensions of the dataset.\u001b[39m\u001b[33m\"\u001b[39m,\n\u001b[32m 395\u001b[39m actual=\u001b[38;5;28mlen\u001b[39m(dims),\n\u001b[32m 396\u001b[39m expected=x.ndim,\n\u001b[32m 397\u001b[39m )\n\u001b[32m 399\u001b[39m new_dims: Sequence[\u001b[38;5;28mstr\u001b[39m] | Sequence[\u001b[38;5;28;01mNone\u001b[39;00m] | \u001b[38;5;28;01mNone\u001b[39;00m\n\u001b[32m 400\u001b[39m \u001b[38;5;28;01mif\u001b[39;00m infer_dims_and_coords:\n",
187-
"\u001b[31mShapeError\u001b[39m: Length of `dims` must match the dimensions of the dataset. (actual 1 != expected 2)"
225+
"\u001b[31mAttributeError\u001b[39m Traceback (most recent call last)",
226+
"\u001b[36mCell\u001b[39m\u001b[36m \u001b[39m\u001b[32mIn[4]\u001b[39m\u001b[32m, line 1\u001b[39m\n\u001b[32m----> \u001b[39m\u001b[32m1\u001b[39m result = \u001b[43mcp\u001b[49m\u001b[43m.\u001b[49m\u001b[43mStructuralTimeSeries\u001b[49m\u001b[43m(\u001b[49m\n\u001b[32m 2\u001b[39m \u001b[43m \u001b[49m\u001b[43mdf\u001b[49m\u001b[43m,\u001b[49m\n\u001b[32m 3\u001b[39m \u001b[43m \u001b[49m\u001b[43mtreatment_time\u001b[49m\u001b[43m,\u001b[49m\n\u001b[32m 4\u001b[39m \u001b[43m \u001b[49m\u001b[43mformula\u001b[49m\u001b[43m=\u001b[49m\u001b[33;43m\"\u001b[39;49m\u001b[33;43my ~ 1 + t + C(month)\u001b[39;49m\u001b[33;43m\"\u001b[39;49m\u001b[43m,\u001b[49m\n\u001b[32m 5\u001b[39m \u001b[43m \u001b[49m\u001b[43mmodel\u001b[49m\u001b[43m=\u001b[49m\u001b[43mcp\u001b[49m\u001b[43m.\u001b[49m\u001b[43mpymc_models\u001b[49m\u001b[43m.\u001b[49m\u001b[43mLinearRegression\u001b[49m\u001b[43m(\u001b[49m\u001b[43msample_kwargs\u001b[49m\u001b[43m=\u001b[49m\u001b[43m{\u001b[49m\u001b[33;43m\"\u001b[39;49m\u001b[33;43mrandom_seed\u001b[39;49m\u001b[33;43m\"\u001b[39;49m\u001b[43m:\u001b[49m\u001b[43m \u001b[49m\u001b[43mseed\u001b[49m\u001b[43m}\u001b[49m\u001b[43m)\u001b[49m\u001b[43m,\u001b[49m\n\u001b[32m 6\u001b[39m \u001b[43m)\u001b[49m\n",
227+
"\u001b[36mFile \u001b[39m\u001b[32m~/Documents/GitHub/CausalPy/causalpy/experiments/structural_time_series.py:145\u001b[39m, in \u001b[36mStructuralTimeSeries.__init__\u001b[39m\u001b[34m(self, data, treatment_time, formula, model, **kwargs)\u001b[39m\n\u001b[32m 138\u001b[39m \u001b[38;5;28;01mif\u001b[39;00m \u001b[38;5;28misinstance\u001b[39m(\u001b[38;5;28mself\u001b[39m.model, PyMCModel):\n\u001b[32m 139\u001b[39m \u001b[38;5;66;03m# Coords for pre-intervention period (for scoring)\u001b[39;00m\n\u001b[32m 140\u001b[39m pre_score_coords = {\n\u001b[32m 141\u001b[39m \u001b[33m\"\u001b[39m\u001b[33mcoeffs\u001b[39m\u001b[33m\"\u001b[39m: \u001b[38;5;28mself\u001b[39m.labels, \u001b[38;5;66;03m# Needed if model uses them for identifying X structure\u001b[39;00m\n\u001b[32m 142\u001b[39m \u001b[33m\"\u001b[39m\u001b[33mdatetime_index\u001b[39m\u001b[33m\"\u001b[39m: \u001b[38;5;28mself\u001b[39m.datapre.index, \u001b[38;5;66;03m# For BSTS\u001b[39;00m\n\u001b[32m 143\u001b[39m \u001b[38;5;66;03m# obs_ind might not be needed by score/predict if model infers from X shape\u001b[39;00m\n\u001b[32m 144\u001b[39m }\n\u001b[32m--> \u001b[39m\u001b[32m145\u001b[39m \u001b[38;5;28mself\u001b[39m.score = \u001b[38;5;28;43mself\u001b[39;49m\u001b[43m.\u001b[49m\u001b[43mmodel\u001b[49m\u001b[43m.\u001b[49m\u001b[43mscore\u001b[49m\u001b[43m(\u001b[49m\n\u001b[32m 146\u001b[39m \u001b[43m \u001b[49m\u001b[43mX\u001b[49m\u001b[43m=\u001b[49m\u001b[43mX_score\u001b[49m\u001b[43m,\u001b[49m\u001b[43m \u001b[49m\u001b[43my\u001b[49m\u001b[43m=\u001b[49m\u001b[38;5;28;43mself\u001b[39;49m\u001b[43m.\u001b[49m\u001b[43mpre_y\u001b[49m\u001b[43m,\u001b[49m\u001b[43m \u001b[49m\u001b[43mcoords\u001b[49m\u001b[43m=\u001b[49m\u001b[43mpre_score_coords\u001b[49m\n\u001b[32m 147\u001b[39m \u001b[43m \u001b[49m\u001b[43m)\u001b[49m\n\u001b[32m 148\u001b[39m \u001b[38;5;28;01melif\u001b[39;00m \u001b[38;5;28misinstance\u001b[39m(\u001b[38;5;28mself\u001b[39m.model, RegressorMixin):\n\u001b[32m 149\u001b[39m \u001b[38;5;28mself\u001b[39m.score = \u001b[38;5;28mself\u001b[39m.model.score(X=X_score, y=\u001b[38;5;28mself\u001b[39m.pre_y)\n",
228+
"\u001b[36mFile \u001b[39m\u001b[32m~/Documents/GitHub/CausalPy/causalpy/pymc_models.py:218\u001b[39m, in \u001b[36mPyMCModel.score\u001b[39m\u001b[34m(self, X, y, coords, **kwargs)\u001b[39m\n\u001b[32m 216\u001b[39m \u001b[38;5;28;01mfor\u001b[39;00m i, unit \u001b[38;5;129;01min\u001b[39;00m \u001b[38;5;28menumerate\u001b[39m(mu_data.coords[\u001b[33m\"\u001b[39m\u001b[33mtreated_units\u001b[39m\u001b[33m\"\u001b[39m].values):\n\u001b[32m 217\u001b[39m unit_mu = mu_data.sel(treated_units=unit).T \u001b[38;5;66;03m# (sample, obs_ind)\u001b[39;00m\n\u001b[32m--> \u001b[39m\u001b[32m218\u001b[39m unit_y = \u001b[43my\u001b[49m\u001b[43m.\u001b[49m\u001b[43msel\u001b[49m(treated_units=unit).data\n\u001b[32m 219\u001b[39m unit_score = r2_score(unit_y, unit_mu.data)\n\u001b[32m 220\u001b[39m scores[\u001b[33mf\u001b[39m\u001b[33m\"\u001b[39m\u001b[33munit_\u001b[39m\u001b[38;5;132;01m{\u001b[39;00mi\u001b[38;5;132;01m}\u001b[39;00m\u001b[33m_r2\u001b[39m\u001b[33m\"\u001b[39m] = unit_score[\u001b[33m\"\u001b[39m\u001b[33mr2\u001b[39m\u001b[33m\"\u001b[39m]\n",
229+
"\u001b[31mAttributeError\u001b[39m: 'numpy.ndarray' object has no attribute 'sel'"
188230
]
189231
}
190232
],

0 commit comments

Comments
 (0)