<jupyter_start><jupyter_text>This is the most basic example from pymc3's "Get started with PyMC3" pageAssume you have a variable mu that is distributed as a normal distrbution,Y ~ N(mu, var) where "~" means is distributed as, mu is the expected value, and var is the variance error of the instrument (std or sigma is the standard deviation). Then assume that mu is a linear function of dependent parameters: alpha, beta1, and beta2, and independent parameters X1 and X2, i.e.                    mu = alpha + beta1*X1 + beta2*X2For this model, assume a weakly informed (uninformed) prior.<jupyter_code>import numpy as np
import matplotlib.pyplot as plt
import pymc3 as pm
from scipy import optimize
%matplotlib inline<jupyter_output>WARNING (theano.tensor.blas): Using NumPy C-API based implementation for BLAS functions.
/packages2/.packages2/x86_64-pc-linux-gnu-rhel7/python/3.6.0/lib/python3.6/site-packages/h5py/__init__.py:36: FutureWarning: Conversion of the second argument of issubdtype from `float` to `np.floating` is deprecated. In future, it will be treated as `np.float64 == np.dtype(float).type`.
  from ._conv import register_converters as _register_converters<jupyter_text>Generate the "true" data<jupyter_code># Intialize random number generator so we get the same set of values after every iteration.
np.random.seed(123) 
# The true dependent parameter values:
alpha = 1; sigma =1; beta = [1, 2.5]
n = 100 # Number of data points

X1 = np.random.randn(n) # Numpy's random normal distribution.
X2 = np.random.randn(n) * 0.2 # scale param
# Simulate observations
Y = alpha + beta[0]*X1 + beta[1]*X2 + np.random.randn(n)*sigma
fig, axes = plt.subplots(1, 2, sharex=True, figsize=(10,4))
axes[0].scatter(X1, Y)
axes[1].scatter(X2, Y)
axes[0].set_ylabel('Y (observation)'); axes[0].set_xlabel('X1 (independent var)'); axes[1].set_xlabel('X2 (independent var)');<jupyter_output><empty_output><jupyter_text>Build the model<jupyter_code>basic_model = pm.Model()

with basic_model: # This is a context manager (all command inside the width statement will be added to pm.Model(), and stay with it)
    # Prior distributions of unknowns for mu. These are STOCASTIC variables.
    alpha = pm.Normal('alpha', mu=0, sd=10)
    beta  = pm.Normal('beta', mu=0, sd=10, shape=2)
    sigma = pm.HalfNormal('sigma', sd=1)
    
    # Define the expected value of outcome (the functional shape of what we will observe, DO NOT confuse this with the 
    # normal distribution that is assosiated with the uncertainty!). 
    # This value is DETERMANISTIC (depends on the parent values that can be stocastic)
    mu = alpha + beta[0]*X1 + beta[1]*X2
    
    # Likelihood (sampling distribution) of observations
    # Observed stocastic, the data likeleyhood, with the observed values of Y must be untouched.
    Y_obs = pm.Normal('Y_obs', mu=mu, sd=sigma, observed=Y) #<jupyter_output><empty_output><jupyter_text>Fit the model  Using the maxium a posteriori (MAP) method. Quick and dirty<jupyter_code>map_estimate = pm.find_MAP(model=basic_model)
map_estimate<jupyter_output>logp = -149.58, ||grad|| = 12.242: 100%|██████████| 19/19 [00:00<00:00, 1329.64it/s]<jupyter_text>Using the automated chooser (No-U-Turn Sampler (NUTS))<jupyter_code>with basic_model:
    trace = pm.sample(500) # Draw 500 posterior samples suing NUTS<jupyter_output>Auto-assigning NUTS sampler...
Initializing NUTS using jitter+adapt_diag...
Multiprocess sampling (2 chains in 2 jobs)
NUTS: [sigma_log__, beta, alpha]
100%|██████████| 1000/1000 [00:00<00:00, 1108.09it/s]<jupyter_text>Why did it do 1000 draws?<jupyter_code>trace['beta']
_ = pm.traceplot(trace)
pm.summary(trace)<jupyter_output><empty_output>
{'fraction_non_alphanumeric': 0.11422355180853958, 'fraction_numerical': 0.03834647810715257, 'mean_word_length': 4.990228013029316, 'pattern_counts': {'":': 0, '<': 24, '<?xml version=': 0, '>': 22, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 0, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'chain_length': '6', 'dir': 'jupyter-structured-clean-dedup', 'id': '127882', 'license': 'BSD-3-Clause', 'n_tokens_mistral': 1341, 'n_tokens_neox': 1232, 'n_words': 397, 'path': 'tutorial_examples/0_linear_regression.ipynb', 'repo_name': 'balarsen/pymc_learning', 'size': '3677'}