0
# Results Analysis
1

2
Tools for accessing, analyzing, and transforming MCMC samples into usable formats. The Fit class provides a dictionary-like interface for working with posterior samples and includes utilities for data export and analysis.
3

4
## Capabilities
5

6
### Sample Access
7

8
Access posterior samples through a dictionary-like interface.
9

10
```python { .api }
11
class Fit:
12
    """
13
    Stores draws from one or more chains.
14
    
15
    A Fit instance works like a Python dictionary. A user-friendly view of draws
16
    is available via to_frame().
17
    
18
    Attributes:
19
        stan_outputs: Raw Stan output for each chain
20
        num_chains: Number of chains
21
        param_names: Parameter names from model
22
        constrained_param_names: All constrained parameter names
23
        dims: Parameter dimensions
24
        num_warmup: Number of warmup iterations
25
        num_samples: Number of sampling iterations
26
        num_thin: Thinning interval
27
        save_warmup: Whether warmup samples were saved
28
        sample_and_sampler_param_names: Names of sample and sampler parameters
29
    """
30
    
31
    def __getitem__(self, param: str):
32
        """
33
        Access parameter draws by name.
34
        
35
        Args:
36
            param: Parameter name
37
            
38
        Returns:
39
            numpy.ndarray: Parameter draws with shape (num_draws, num_chains)
40
            
41
        Raises:
42
            KeyError: If parameter name not found
43
        """
44
    
45
    def __contains__(self, key: str) -> bool:
46
        """
47
        Check if parameter exists in the fit.
48
        
49
        Args:
50
            key: Parameter name
51
            
52
        Returns:
53
            bool: True if parameter exists
54
        """
55
    
56
    def __iter__(self):
57
        """
58
        Iterate over parameter names.
59
        
60
        Yields:
61
            str: Parameter names in order
62
        """
63
    
64
    def __len__(self) -> int:
65
        """
66
        Number of parameters in the fit.
67
        
68
        Returns:
69
            int: Total number of parameters
70
        """
71
```
72

73
### Data Export
74

75
Convert samples to user-friendly formats.
76

77
```python { .api }
78
def to_frame(self):
79
    """
80
    Return view of draws as a pandas DataFrame.
81
    
82
    The DataFrame contains all parameters and diagnostic information,
83
    flattened across all chains with draws as rows and parameters as columns.
84
    
85
    Returns:
86
        pandas.DataFrame: DataFrame with num_draws rows and 
87
                         num_flat_params columns
88
                         
89
    Raises:
90
        RuntimeError: If pandas is not installed
91
        
92
    Notes:
93
        - Requires pandas to be installed
94
        - Includes both model parameters and sampler diagnostics
95
        - All chains are combined into a single DataFrame
96
        - Column names match sample_and_sampler_param_names + constrained_param_names
97
    """
98
```
99

100
### String Representation
101

102
Human-readable summary of fit contents.
103

104
```python { .api }
105
def __repr__(self) -> str:
106
    """
107
    String representation of the Fit showing parameter summaries.
108
    
109
    Provides a concise overview of all parameters including:
110
    - Parameter names and dimensions
111
    - Number of chains and draws
112
    - Brief statistical summary
113
    
114
    Returns:
115
        str: Formatted summary of fit contents
116
    """
117
```
118

119
## Usage Examples
120

121
### Basic Sample Access
122

123
```python
124
import stan
125

126
program_code = """
127
parameters {
128
    real mu;
129
    real<lower=0> sigma;
130
    vector[3] theta;
131
}
132
model {
133
    mu ~ normal(0, 1);
134
    sigma ~ exponential(1);
135
    theta ~ normal(0, 1);
136
}
137
"""
138

139
model = stan.build(program_code)
140
fit = model.sample(num_chains=4, num_samples=1000)
141

142
# Access individual parameters
143
mu_samples = fit['mu']
144
sigma_samples = fit['sigma']
145
theta_samples = fit['theta']
146

147
print(f"mu samples shape: {mu_samples.shape}")
148
print(f"sigma samples shape: {sigma_samples.shape}")
149
print(f"theta samples shape: {theta_samples.shape}")
150

151
# Check parameter existence
152
print(f"'mu' in fit: {'mu' in fit}")
153
print(f"'nonexistent' in fit: {'nonexistent' in fit}")
154

155
# Iterate over parameters
156
print("All parameters:")
157
for param_name in fit:
158
    print(f"  {param_name}: {fit[param_name].shape}")
159

160
print(f"Total parameters: {len(fit)}")
161
```
162

163
### DataFrame Export and Analysis
164

165
```python
166
import stan
167
import numpy as np
168
import pandas as pd
169

170
program_code = """
171
data {
172
    int<lower=0> N;
173
    vector[N] x;
174
    vector[N] y;
175
}
176
parameters {
177
    real alpha;
178
    real beta;
179
    real<lower=0> sigma;
180
}
181
model {
182
    alpha ~ normal(0, 10);
183
    beta ~ normal(0, 10);
184
    sigma ~ exponential(1);
185
    y ~ normal(alpha + beta * x, sigma);
186
}
187
"""
188

189
# Generate data
190
N = 50
191
x = np.random.normal(0, 1, N)
192
y = 2 + 3 * x + np.random.normal(0, 1, N)
193
data = {'N': N, 'x': x.tolist(), 'y': y.tolist()}
194

195
model = stan.build(program_code, data=data)
196
fit = model.sample(num_chains=4, num_samples=1000)
197

198
# Convert to DataFrame
199
df = fit.to_frame()
200
print(f"DataFrame shape: {df.shape}")
201
print(f"DataFrame columns: {list(df.columns)}")
202

203
# Basic statistics
204
print("\nParameter summaries:")
205
print(df.describe())
206

207
# Chain-specific analysis
208
print("\nMean by chain:")
209
print(df.groupby(level=1, axis=1).mean())
210

211
# Plot traces (requires matplotlib)
212
import matplotlib.pyplot as plt
213

214
fig, axes = plt.subplots(3, 1, figsize=(10, 8))
215
parameters = ['alpha', 'beta', 'sigma']
216

217
for i, param in enumerate(parameters):
218
    for chain in range(4):
219
        axes[i].plot(df[(param, chain)], alpha=0.7, label=f'Chain {chain}')
220
    axes[i].set_title(f'{param} trace')
221
    axes[i].legend()
222

223
plt.tight_layout()
224
plt.show()
225
```
226

227
### Diagnostic Analysis
228

229
```python
230
import stan
231
import numpy as np
232

233
program_code = """
234
parameters {
235
    real mu;
236
    real<lower=0> sigma;
237
}
238
model {
239
    mu ~ normal(0, 1);
240
    sigma ~ exponential(1);
241
}
242
"""
243

244
model = stan.build(program_code)
245
fit = model.sample(num_chains=4, num_samples=1000, num_warmup=1000)
246

247
# Access diagnostic parameters
248
print("Available parameters:")
249
for param in fit:
250
    print(f"  {param}")
251

252
# Check sampler diagnostics
253
if 'accept_stat__' in fit:
254
    accept_stat = fit['accept_stat__']
255
    print(f"\nAcceptance rate statistics:")
256
    print(f"  Mean: {np.mean(accept_stat):.3f}")
257
    print(f"  Min: {np.min(accept_stat):.3f}")
258
    print(f"  Max: {np.max(accept_stat):.3f}")
259

260
if 'treedepth__' in fit:
261
    treedepth = fit['treedepth__']
262
    print(f"\nTree depth statistics:")
263
    print(f"  Mean: {np.mean(treedepth):.1f}")
264
    print(f"  Max: {np.max(treedepth)}")
265

266
if 'stepsize__' in fit:
267
    stepsize = fit['stepsize__']
268
    print(f"\nStep size by chain:")
269
    for chain in range(fit.num_chains):
270
        chain_stepsize = stepsize[:, chain]
271
        print(f"  Chain {chain}: {np.mean(chain_stepsize):.4f}")
272

273
# Print fit summary
274
print(f"\nFit summary:")
275
print(fit)
276
```
277

278
### Advanced Analysis
279

280
```python
281
import stan
282
import numpy as np
283
from scipy import stats
284

285
program_code = """
286
parameters {
287
    real mu;
288
    real<lower=0> sigma;
289
}
290
model {
291
    mu ~ normal(0, 1);
292
    sigma ~ exponential(1);
293
}
294
"""
295

296
model = stan.build(program_code)
297
fit = model.sample(num_chains=4, num_samples=2000)
298

299
# Posterior analysis
300
mu_samples = fit['mu'].flatten()  # Combine all chains
301
sigma_samples = fit['sigma'].flatten()
302

303
print("Posterior summaries:")
304
print(f"mu: mean={np.mean(mu_samples):.3f}, std={np.std(mu_samples):.3f}")
305
print(f"sigma: mean={np.mean(sigma_samples):.3f}, std={np.std(sigma_samples):.3f}")
306

307
# Credible intervals
308
mu_ci = np.percentile(mu_samples, [2.5, 97.5])
309
sigma_ci = np.percentile(sigma_samples, [2.5, 97.5])
310

311
print(f"\n95% Credible Intervals:")
312
print(f"mu: [{mu_ci[0]:.3f}, {mu_ci[1]:.3f}]")
313
print(f"sigma: [{sigma_ci[0]:.3f}, {sigma_ci[1]:.3f}]")
314

315
# Chain convergence (R-hat approximation)
316
def split_rhat(chains):
317
    """Simple R-hat calculation for demonstration"""
318
    n_chains, n_draws = chains.shape[1], chains.shape[0]
319
    
320
    # Split each chain in half
321
    first_half = chains[:n_draws//2, :]
322
    second_half = chains[n_draws//2:, :]
323
    
324
    # Combine split chains
325
    all_chains = np.concatenate([first_half, second_half], axis=1)
326
    
327
    # Between and within chain variance
328
    chain_means = np.mean(all_chains, axis=0)
329
    overall_mean = np.mean(all_chains)
330
    
331
    B = (n_draws // 2) * np.var(chain_means, ddof=1)
332
    W = np.mean([np.var(all_chains[:, i], ddof=1) for i in range(all_chains.shape[1])])
333
    
334
    var_est = ((n_draws // 2 - 1) / (n_draws // 2)) * W + B / (n_draws // 2)
335
    rhat = np.sqrt(var_est / W)
336
    
337
    return rhat
338

339
mu_rhat = split_rhat(fit['mu'])
340
sigma_rhat = split_rhat(fit['sigma'])
341

342
print(f"\nR-hat diagnostics:")
343
print(f"mu: {mu_rhat:.3f}")
344
print(f"sigma: {sigma_rhat:.3f}")
345
print("(Values < 1.1 indicate good convergence)")
346
```