### Install pymbar with Conda or Pip

Source: https://context7.com/choderalab/pymbar/llms.txt

Instructions for installing the pymbar library using Conda or Pip. It provides options for installing with JAX acceleration (recommended) or without JAX for a smaller footprint.

```bash
conda install -c conda-forge pymbar

pip install pymbar[jax]

pip install pymbar
```

--------------------------------

### Default Solver Protocol Example in Python

Source: https://github.com/choderalab/pymbar/blob/main/docs/strategies_for_solution.rst

Demonstrates the construction of a solver protocol using a tuple of dictionaries, where each dictionary represents an optimization operation. This example shows the 'default' protocol, which aims for high accuracy by using the 'hybr' method with continuation.

```python
solver_protocol = (
        dict(method="hybr", continuation=True),
   )
```

--------------------------------

### Install PyMBAR from GitHub

Source: https://github.com/choderalab/pymbar/blob/main/README.md

Installs the PyMBAR package directly from its GitHub repository using pip. This is useful for obtaining the latest development version.

```bash
pip install git+https://github.com/choderalab/pymbar.git
```

--------------------------------

### Disable JAX Acceleration with Pip

Source: https://github.com/choderalab/pymbar/blob/main/README.md

Installs the `pymbar` package using pip. In versions >= 4.2.0, this installation method respects the `PYMBAR_DISABLE_JAX` environment variable.

```bash
pip install pymbar
```

--------------------------------

### Disable JAX Acceleration with Conda

Source: https://github.com/choderalab/pymbar/blob/main/README.md

Installs the `pymbar-core` package from conda-forge, which does not include JAX as a dependency, effectively disabling JAX acceleration.

```bash
conda install -c conda-forge pymbar-core
```

--------------------------------

### Disable JAX Acceleration via Environment Variable

Source: https://github.com/choderalab/pymbar/blob/main/README.md

Sets the `PYMBAR_DISABLE_JAX` environment variable to a truthy value to disable JAX acceleration, even if the `jax` package is installed. This method is available in `pymbar` versions 4.2.0 and later.

```bash
export PYMBAR_DISABLE_JAX=TRUE
```

--------------------------------

### Detect Equilibration and Subsample Correlated Data (Python)

Source: https://github.com/choderalab/pymbar/blob/main/docs/timeseries.rst

Detects the equilibration point in a timeseries and subsamples it to obtain uncorrelated data. It returns the start index of the production region, the statistical inefficiency, and the maximum number of uncorrelated samples. The production region is then extracted and subsampled.

```python
from pymbar import timeseries
t0, g, Neff_max = timeseries.detect_equilibration(A_t) # compute indices of uncorrelated timeseries
A_t_equil = A_t[t0:]
indices = timeseries.subsample_correlated_data(A_t_equil, g=g)
A_n = A_t_equil[indices]
```

--------------------------------

### Sample Harmonic Oscillators and Initialize MBAR (Python)

Source: https://github.com/choderalab/pymbar/blob/main/README.md

Demonstrates sampling data from a 1D harmonic oscillator test system and initializing the MBAR object. It requires the `pymbar` library and its test systems module.

```python
from pymbar import testsystems
x_n, u_kn, N_k, s_n = testsystems.HarmonicOscillatorsTestCase().sample()
```

```python
from pymbar import MBAR
mbar = MBAR(u_kn, N_k)
```

--------------------------------

### Initialize MBAR and Access Free Energies and Weights (Python)

Source: https://context7.com/choderalab/pymbar/llms.txt

Demonstrates how to initialize the MBAR class with reduced potential energies and sample counts, and how to access the computed dimensionless free energies and the weight matrix.

```python
import numpy as np
from pymbar import MBAR, testsystems

# Generate test data from harmonic oscillators at different temperatures
x_n, u_kn, N_k, s_n = testsystems.HarmonicOscillatorsTestCase().sample(mode='u_kn')

# Initialize MBAR with reduced potential matrix and sample counts
# u_kn[k,n] is the reduced potential of sample n evaluated at state k
# N_k[k] is the number of samples from state k
mbar = MBAR(u_kn, N_k)

# Access computed dimensionless free energies
print(f"Dimensionless free energies: {mbar.f_k}")

# Access weight matrix (N samples x K states)
weights = mbar.W_nk  # or mbar.weights()
print(f"Weight matrix shape: {weights.shape}")
```

--------------------------------

### Compute Free Energy Differences and Uncertainties (Python)

Source: https://context7.com/choderalab/pymbar/llms.txt

Shows how to compute the matrix of free energy differences between all pairs of thermodynamic states, including their statistical uncertainties. It also demonstrates how to obtain the covariance matrix for log weights.

```python
from pymbar import MBAR, testsystems

# Generate and analyze harmonic oscillator data
x_n, u_kn, N_k, s_n = testsystems.HarmonicOscillatorsTestCase().sample(mode='u_kn')
mbar = MBAR(u_kn, N_k)

# Compute free energy differences with uncertainties
results = mbar.compute_free_energy_differences()

Delta_f = results['Delta_f']   # K x K matrix of free energy differences
dDelta_f = results['dDelta_f'] # K x K matrix of uncertainties

# Free energy difference between states 0 and 4
print(f"Delta F(0->4) = {Delta_f[0,4]:.4f} +/- {dDelta_f[0,4]:.4f} kT")

# Free energy differences are antisymmetric: Delta_f[i,j] = -Delta_f[j,i]
print(f"Delta F(4->0) = {Delta_f[4,0]:.4f} kT")

# Request the covariance matrix for error propagation
results_with_theta = mbar.compute_free_energy_differences(return_theta=True)
Theta = results_with_theta['Theta']  # Covariance matrix for log weights
```

--------------------------------

### Compute Free Energies and Uncertainties with pymbar

Source: https://context7.com/choderalab/pymbar/llms.txt

This snippet demonstrates how to compute free energy differences using pymbar, comparing bootstrap and analytical uncertainty estimation methods. It also shows how to compute expectation values and perform entropy/enthalpy decomposition with bootstrap uncertainties.

```python
from pymbar import MBAR

# Assume u_kn, N_k, x_n are defined elsewhere
# u_kn: Reduced potential energies
# N_k: Number of samples in each state
# x_n: Observables

# Compute free energies with bootstrap uncertainties
results = mbar.compute_free_energy_differences(uncertainty_method='bootstrap')
Delta_f = results['Delta_f']
dDelta_f_boot = results['dDelta_f']

# Compare with analytical uncertainties
results_analytical = mbar.compute_free_energy_differences(uncertainty_method='svd-ew')
dDelta_f_anal = results_analytical['dDelta_f']

print("Uncertainty comparison (states 0->4):")
print(f"  Bootstrap:  {dDelta_f_boot[0,4]:.4f}")
print(f"  Analytical: {dDelta_f_anal[0,4]:.4f}")

# Bootstrap expectations
results_exp = mbar.compute_expectations(x_n, uncertainty_method='bootstrap')
print(f"\n<x> at state 0: {results_exp['mu'][0]:.4f} +/- {results_exp['sigma'][0]:.4f}")

# Bootstrap entropy/enthalpy decomposition
results_ent = mbar.compute_entropy_and_enthalpy(uncertainty_method='bootstrap')
```

--------------------------------

### Compute Perturbed Free Energies for New States (Python)

Source: https://context7.com/choderalab/pymbar/llms.txt

Calculates free energies at new thermodynamic states not present in the original sampling. This allows extrapolation using existing MBAR weights. Requires an MBAR object and the reduced potentials for the new states.

```python
from pymbar import MBAR, testsystems
import numpy as np

# Generate data from original states
x_n, u_kn, N_k, s_n = testsystems.HarmonicOscillatorsTestCase().sample(mode='u_kn')
mbar = MBAR(u_kn, N_k)

# Define new states (must evaluate samples at new potentials)
# Here we use the original states as an example
u_ln = u_kn  # L x N matrix of new reduced potentials

results = mbar.compute_perturbed_free_energies(u_ln)

Delta_f = results['Delta_f']   # L x L matrix of free energy differences
dDelta_f = results['dDelta_f'] # L x L matrix of uncertainties

print(f"Free energy difference (new states 0->4): {Delta_f[0,4]:.4f} +/- {dDelta_f[0,4]:.4f}")

# Example: Compute free energy at an intermediate temperature
# (assuming linear interpolation of potential)
alpha = 0.5
u_new = alpha * u_kn[2, :] + (1 - alpha) * u_kn[3, :]
u_ln_new = u_new.reshape(1, -1)  # Single new state
results_new = mbar.compute_perturbed_free_energies(
    np.vstack([u_kn, u_ln_new])  # Append new state to existing
)
```

--------------------------------

### Compute Observable Averages and Differences (Python)

Source: https://context7.com/choderalab/pymbar/llms.txt

Illustrates how to compute expectation values of observables at each thermodynamic state using MBAR reweighting. It also covers computing differences in observables between states and handling state-dependent observables.

```python
from pymbar import MBAR, testsystems
import numpy as np

# Generate data
x_n, u_kn, N_k, s_n = testsystems.HarmonicOscillatorsTestCase().sample(mode='u_kn')
mbar = MBAR(u_kn, N_k)

# Compute expectation of position at each state
A_n = x_n  # Observable: position
results = mbar.compute_expectations(A_n)

mu = results['mu']      # Expected values at each state
sigma = results['sigma'] # Standard errors

for k in range(len(mu)):
    print(f"State {k}: <x> = {mu[k]:.4f} +/- {sigma[k]:.4f}")

# Compute differences in observable between states
results_diff = mbar.compute_expectations(A_n, output='differences')
Delta_A = results_diff['mu']   # K x K matrix of observable differences
dDelta_A = results_diff['sigma']

print(f"<x>_4 - <x>_0 = {Delta_A[0,4]:.4f} +/- {dDelta_A[0,4]:.4f}")

# State-dependent observables (different observable at each state)
A_kn = u_kn  # Use energies as state-dependent observable
results_sd = mbar.compute_expectations(A_kn, state_dependent=True)
```

--------------------------------

### Construct Free Energy Surfaces (FES) with MBAR

Source: https://context7.com/choderalab/pymbar/llms.txt

Computes potentials of mean force from biased simulation data. Supports histogram-based construction and bootstrap uncertainty estimation.

```python
from pymbar import FES, testsystems
import numpy as np

# Generate umbrella sampling data
x_n, u_kn, N_k, s_n = testsystems.HarmonicOscillatorsTestCase().sample(mode='u_kn', seed=0)

# Initialize FES with MBAR
fes = FES(u_kn, N_k)
u_n = u_kn[0, :]

# Generate FES with bootstrap uncertainty
N_tot = N_k.sum()
x_sorted = np.sort(x_n)
bins = np.append(x_sorted[0::int(N_tot/15)], x_sorted.max() + 0.1)
histogram_params = {'bin_edges': [bins]}

fes.generate_fes(u_n, x_n, fes_type='histogram',
                 histogram_parameters=histogram_params,
                 n_bootstraps=50, seed=42)

# Get FES with bootstrap uncertainties
x_query = np.linspace(x_n.min() + 0.1, x_n.max() - 0.1, 30)
results = fes.get_fes(x_query, uncertainty_method='bootstrap')

f_i = results['f_i']
df_i = results['df_i']
```

--------------------------------

### Compute Expectation Values with PyMBAR (Python)

Source: https://github.com/choderalab/pymbar/blob/main/README.md

Estimates expectation values and their uncertainties for observables across different thermodynamic states using an MBAR object. Requires observable data in a specific format.

```python
A_kn = x_kn # use position of harmonic oscillator as observable
results = mbar.compute_expectations(A_kn)
```

--------------------------------

### Enabling Logging in pymbar

Source: https://github.com/choderalab/pymbar/blob/main/docs/moving_from_pymbar3.rst

Demonstrates how to configure Python's logging module to capture output from pymbar. This replaces the previous behavior of printing directly to standard output when verbose mode was enabled, offering more control over message handling.

```python
import logging
import sys
logging.basicConfig(stream=sys.stdout, level=logging.INFO)
```

--------------------------------

### Configure Bootstrap Solver Protocol (Python)

Source: https://github.com/choderalab/pymbar/blob/main/docs/strategies_for_solution.rst

Defines the configuration for the bootstrap solver protocol using an adaptive method. It sets a tolerance of 1e-6 and specifies options for minimum sub-chain iterations and maximum iterations. This configuration is suitable for calculating variances from bootstrapped data where high precision is not always required.

```python
bootstrap_solver_protocol = (dict(method="adaptive", tol = 1e-6, options=dict(min_sc_iter=0,maximum_iterations=100)))
```

--------------------------------

### Compute Entropy and Enthalpy Decomposition (Python)

Source: https://context7.com/choderalab/pymbar/llms.txt

Decomposes free energy differences into entropic and enthalpic contributions, providing insight into the driving forces of thermodynamic changes. Requires an MBAR object. Outputs free energy, enthalpy, and entropy differences with their uncertainties.

```python
from pymbar import MBAR, testsystems

# Generate data
x_n, u_kn, N_k, s_n = testsystems.HarmonicOscillatorsTestCase().sample(mode='u_kn')
mbar = MBAR(u_kn, N_k)

# Compute entropy and enthalpy decomposition
results = mbar.compute_entropy_and_enthalpy()

Delta_f = results['Delta_f']   # Free energy differences (kT)
Delta_u = results['Delta_u']   # Enthalpy/energy differences (kT)
Delta_s = results['Delta_s']   # Entropy differences (S/k_B)

dDelta_f = results['dDelta_f'] # Uncertainties
dDelta_u = results['dDelta_u']
dDelta_s = results['dDelta_s']

# Thermodynamic identity: Delta_f = Delta_u - Delta_s (at constant T)
i, j = 0, 4
print(f"Delta F = {Delta_f[i,j]:.4f} +/- {dDelta_f[i,j]:.4f} kT")
print(f"Delta U = {Delta_u[i,j]:.4f} +/- {dDelta_u[i,j]:.4f} kT")
print(f"Delta S = {Delta_s[i,j]:.4f} +/- {dDelta_s[i,j]:.4f} k_B")
print(f"Delta U - Delta S = {Delta_u[i,j] - Delta_s[i,j]:.4f} (should equal Delta F)")
```

--------------------------------

### Bootstrap Uncertainty Estimation with MBAR (Python)

Source: https://context7.com/choderalab/pymbar/llms.txt

Demonstrates initializing the MBAR class with bootstrap sampling enabled for uncertainty estimation. This is an alternative to analytical covariance calculations, useful for complex correlation structures or when analytical uncertainties are unreliable. Requires PyMBAR and NumPy.

```python
from pymbar import MBAR, testsystems
import numpy as np

# Generate data
x_n, u_kn, N_k, s_n = testsystems.HarmonicOscillatorsTestCase().sample(mode='u_kn')

# Initialize MBAR with bootstrap sampling
mbar = MBAR(u_kn, N_k, n_bootstraps=200)
```

--------------------------------

### Compute Free Energy Differences with PyMBAR (Python)

Source: https://github.com/choderalab/pymbar/blob/main/README.md

Calculates dimensionless free energy differences and their uncertainties between thermodynamic states using an initialized MBAR object. The results are returned in a dictionary.

```python
results = mbar.compute_free_energy_differences()
```

--------------------------------

### Compute Free Energy with Bennett Acceptance Ratio (BAR)

Source: https://context7.com/choderalab/pymbar/llms.txt

Calculates free energy differences between two states using the BAR method. It supports multiple solver methods and provides tools for assessing phase space overlap.

```python
from pymbar import bar, testsystems

# Generate Gaussian work distributions for forward and reverse processes
w_F, w_R = testsystems.gaussian_work_example(mu_F=None, DeltaF=1.0, seed=42)

# Compute free energy difference using BAR
results = bar(w_F, w_R)

Delta_f = results['Delta_f']
dDelta_f = results['dDelta_f']

print(f"Free energy difference: {Delta_f:.4f} +/- {dDelta_f:.4f} kT")

# Alternative methods for solving BAR equation
results_bisect = bar(w_F, w_R, method='bisection')
results_sci = bar(w_F, w_R, method='self-consistent-iteration')

# Check overlap between forward and reverse distributions
from pymbar import bar_overlap
overlap = bar_overlap(w_F, w_R)
print(f"Phase space overlap: {overlap:.4f}")
```

--------------------------------

### Detect Equilibration Time in Timeseries (Python)

Source: https://context7.com/choderalab/pymbar/llms.txt

Automatically detects the equilibration time of a simulation timeseries by maximizing the number of effective uncorrelated samples in the production phase. Offers both standard and binary search methods for detection. Requires NumPy.

```python
from pymbar import timeseries
from pymbar.testsystems import correlated_timeseries_example
import numpy as np

# Create timeseries with initial equilibration period
# Simulate: 1000 equilibration steps + 9000 production steps
equilibration = np.linspace(5.0, 0.0, 1000)  # Drift during equilibration
production = correlated_timeseries_example(N=9000, tau=5.0)
A_t = np.concatenate([equilibration, production])

# Detect equilibration automatically
t0, g, Neff = timeseries.detect_equilibration(A_t)

print(f"Equilibration ends at sample: {t0}")
print(f"Statistical inefficiency of production: {g:.2f}")
print(f"Effective uncorrelated samples: {Neff:.0f}")
# Output:
# Equilibration ends at sample: ~1000
# Statistical inefficiency of production: ~11
# Effective uncorrelated samples: ~818

# Extract production data
A_production = A_t[t0:]
print(f"Production samples: {len(A_production)}")

# Use binary search for faster detection on large datasets
t0_bs, g_bs, Neff_bs = timeseries.detect_equilibration_binary_search(A_t)
print(f"Binary search result: t0={t0_bs}, Neff={Neff_bs:.0f}")

# Subsample production data
indices = timeseries.subsample_correlated_data(A_production, g=g)
A_uncorrelated = A_production[indices]
print(f"Final uncorrelated samples: {len(A_uncorrelated)}")
```

--------------------------------

### Configure Solver Protocol for MBAR

Source: https://github.com/choderalab/pymbar/blob/main/docs/strategies_for_solution.rst

Defines a sequence of optimization methods for solving MBAR equations. It allows specifying different methods like 'hybr', 'L-BFGS-B', and 'adaptive' with their respective tolerances and options. This enables a robust solution process by trying multiple solvers if the initial ones fail.

```python
solver_protocol = (
        dict(method="hybr"),
        dict(method="L-BFGS-B", tol = 1.0e-5, continuation = True, options=dict(maxiter=1000)),
        dict(method="adaptive", tol = 1.0e-12, options=dict(maxiter=1000,min_sc_iter=5))
)
```

--------------------------------

### Compute Overlap Matrix for Sampling Quality (Python)

Source: https://context7.com/choderalab/pymbar/llms.txt

Calculates the overlap matrix between thermodynamic states to assess sampling quality and phase space overlap. Poor overlap indicates unreliable free energy estimates. Requires an MBAR object.

```python
from pymbar import MBAR, testsystems

# Generate data
x_n, u_kn, N_k, s_n = testsystems.HarmonicOscillatorsTestCase().sample(mode='u_kn')
mbar = MBAR(u_kn, N_k)

# Compute overlap metrics
results = mbar.compute_overlap()

overlap_scalar = results['scalar']      # 1 minus second largest eigenvalue
eigenvalues = results['eigenvalues']    # Sorted eigenvalues (descending)
overlap_matrix = results['matrix']      # K x K overlap probability matrix

print(f"Overlap scalar (0=poor, 1=excellent): {overlap_scalar:.4f}")
# Output: Overlap scalar (0=poor, 1=excellent): 0.1142

# Overlap matrix: O[i,j] is probability of observing sample from state i in state j
print("Overlap matrix (first 4x4 block):")
print(overlap_matrix[:4, :4].round(3))
# Each row sums to 1.0

# Check effective sample number at each state
N_eff = mbar.compute_effective_sample_number()
print(f"Effective samples per state: {N_eff.round(1)}")
```

--------------------------------

### Estimate Free Energy with Exponential Averaging

Source: https://context7.com/choderalab/pymbar/llms.txt

Utilizes unidirectional exponential averaging (Zwanzig relation) for quick free energy estimates. Includes support for Gaussian approximations and correlated timeseries data.

```python
from pymbar import exp, exp_gauss, testsystems

# Generate work values
w_F, w_R = testsystems.gaussian_work_example(mu_F=None, DeltaF=1.0, seed=42)

# Forward exponential averaging
results_F = exp(w_F)
print(f"Forward: Delta F = {results_F['Delta_f']:.4f} +/- {results_F['dDelta_f']:.4f} kT")

# Gaussian approximation (faster, assumes Gaussian work distribution)
results_gauss = exp_gauss(w_F)
print(f"Gaussian approx: Delta F = {results_gauss['Delta_f']:.4f} +/- {results_gauss['dDelta_f']:.4f} kT")

# Handle correlated timeseries data
results_ts = exp(w_F, is_timeseries=True)
```

--------------------------------

### Accessing Free Energy Results from EXP Estimator in pymbar

Source: https://github.com/choderalab/pymbar/blob/main/docs/moving_from_pymbar3.rst

Illustrates retrieving free energy difference results from the EXP estimator in pymbar. Similar to BAR, the EXP estimator now returns a dictionary, providing structured access to forward and reverse free energy differences and their uncertainties.

```python
results = exp(w_F)
print(f"Forward free energy difference is {results['Delta_f']:.3f} +- {results['dDelta_f']:.3f} kT)
results = exp(w_R)
print(f"Reverse free energy difference is {results['Delta_f']:.3f} +- {results['dDelta_f']:.3f} kT)
```

--------------------------------

### Compute Multiple Expectations with Covariance (Python)

Source: https://context7.com/choderalab/pymbar/llms.txt

Computes expectations of multiple observables and their covariances simultaneously at a single state. Useful for analyzing correlations between different properties from the same simulation data. Requires MBAR object and potentials at the target state.

```python
from pymbar import MBAR, testsystems
import numpy as np

# Generate data
x_n, u_kn, N_k, s_n = testsystems.HarmonicOscillatorsTestCase().sample(mode='u_kn')
mbar = MBAR(u_kn, N_k)

# Define multiple observables: x, x^2, x^3
A_in = np.array([x_n, x_n**2, x_n**3])

# Compute expectations at state 0
u_n = u_kn[0, :]  # Reduced potential at target state
results = mbar.compute_multiple_expectations(A_in, u_n, compute_covariance=True)

mu = results['mu']           # Expectations: [<x>, <x^2>, <x^3>]
sigma = results['sigma']     # Standard errors
cov = results['covariances'] # Covariance matrix

print(f"<x>   = {mu[0]:.4f} +/- {sigma[0]:.4f}")
print(f"<x^2> = {mu[1]:.4f} +/- {sigma[1]:.4f}")
print(f"<x^3> = {mu[2]:.4f} +/- {sigma[2]:.4f}")
# Output:
# <x>   = 0.0012 +/- 0.0089
# <x^2> = 0.5023 +/- 0.0112
# <x^3> = 0.0041 +/- 0.0267

# Variance can be computed from moments
variance = mu[1] - mu[0]**2
print(f"Var(x) = {variance:.4f}")
```

--------------------------------

### Detect Equilibration with Skipping (Python)

Source: https://github.com/choderalab/pymbar/blob/main/docs/timeseries.rst

Detects the equilibration point in a timeseries, but optimizes the process by only considering every nskip samples as potential time origins. This can speed up computation for long timeseries at the cost of potentially discarding slightly more data.

```python
nskip = 10 # only try every 10 samples for time origins
t0, g, Neff_max = timeseries.detect_equilibration(A_t, nskip=nskip)
```

--------------------------------

### Accessing Free Energy Results from BAR Estimator in pymbar

Source: https://github.com/choderalab/pymbar/blob/main/docs/moving_from_pymbar3.rst

Shows how to retrieve free energy difference results from the BAR estimator in pymbar. The function now returns a dictionary, allowing access to 'Delta_f' and 'dDelta_f' by key, improving code readability and robustness compared to tuple indexing.

```python
results = bar(w_F, w_R)
print(f'Free energy difference is {results['Delta_f']:.3f} +- {results['Delta_f']:.3f} kT')
```

--------------------------------

### Disable JAX Acceleration in pymbar

Source: https://context7.com/choderalab/pymbar/llms.txt

This snippet shows how to disable JAX acceleration in pymbar using an environment variable. It demonstrates how to check if JAX is being used and confirms that MBAR functionality remains consistent.

```python
import os

# Method 1: Environment variable (must be set before importing pymbar)
os.environ['PYMBAR_DISABLE_JAX'] = 'TRUE'

# Now import pymbar
from pymbar import MBAR, testsystems

# Check if JAX is being used
from pymbar import mbar_solvers
print(f"Using JAX: {mbar_solvers.using_jax()}")

# Generate test data
x_n, u_kn, N_k, s_n = testsystems.HarmonicOscillatorsTestCase().sample(mode='u_kn')

# MBAR works the same way regardless of JAX status
mbar = MBAR(u_kn, N_k)
results = mbar.compute_free_energy_differences()
print(f"Delta F(0->4) = {results['Delta_f'][0,4]:.4f} kT")
```

--------------------------------

### Subsample Correlated Data to Extract Uncorrelated Samples (Python)

Source: https://context7.com/choderalab/pymbar/llms.txt

Determines indices for extracting uncorrelated samples from a correlated timeseries, essential for methods like MBAR. Can use pre-computed statistical inefficiency or a conservative uniform spacing. Handles multiple trajectories with a shared statistical inefficiency value. Requires NumPy.

```python
from pymbar import timeseries
from pymbar.testsystems import correlated_timeseries_example
import numpy as np

# Generate correlated timeseries
tau = 15.0
A_t = correlated_timeseries_example(N=50000, tau=tau)
print(f"Original timeseries length: {len(A_t)}")

# Get indices of uncorrelated samples
indices = timeseries.subsample_correlated_data(A_t)
A_n = A_t[indices]  # Uncorrelated subset
print(f"Uncorrelated samples: {len(A_n)}")
# Output: Uncorrelated samples: ~1613

# Use pre-computed statistical inefficiency
g = timeseries.statistical_inefficiency(A_t, fast=True)
indices_g = timeseries.subsample_correlated_data(A_t, g=g)
print(f"Samples with g={g:.1f}: {len(indices_g)}")

# Conservative subsampling (uniform spacing)
indices_cons = timeseries.subsample_correlated_data(A_t, conservative=True)
print(f"Conservative samples: {len(indices_cons)}")

# Process multiple trajectories with same g
trajectories = [correlated_timeseries_example(N=n, tau=tau)
                for n in [10000, 20000, 30000]]
g_combined = timeseries.statistical_inefficiency_multiple(trajectories)
uncorrelated_data = [traj[timeseries.subsample_correlated_data(traj, g=g_combined)]
                     for traj in trajectories]
```

--------------------------------

### Compute Statistical Inefficiency of Timeseries (Python)

Source: https://context7.com/choderalab/pymbar/llms.txt

Calculates the statistical inefficiency of a timeseries, indicating how many correlated samples are equivalent to one uncorrelated sample. Supports standard, fast, and FFT-based methods, as well as cross-correlation between two timeseries. Requires NumPy.

```python
from pymbar import timeseries
from pymbar.testsystems import correlated_timeseries_example
import numpy as np

# Generate correlated timeseries with known correlation time
tau = 10.0  # Correlation time
A_t = correlated_timeseries_example(N=10000, tau=tau)

# Compute statistical inefficiency (g = 1 + 2*tau for exponential correlation)
g = timeseries.statistical_inefficiency(A_t)
print(f"Statistical inefficiency g = {g:.2f} (expected ~{1 + 2*tau:.0f})")
# Output: Statistical inefficiency g = 20.89 (expected ~21)

# Fast mode (less accurate but quicker)
g_fast = timeseries.statistical_inefficiency(A_t, fast=True)
print(f"Fast estimate: g = {g_fast:.2f}")

# FFT-based method (requires statsmodels)
try:
    g_fft = timeseries.statistical_inefficiency(A_t, fft=True)
    print(f"FFT estimate: g = {g_fft:.2f}")
except ImportError:
    print("FFT method requires statsmodels package")

# Cross-correlation between two timeseries
B_t = correlated_timeseries_example(N=10000, tau=tau)
g_cross = timeseries.statistical_inefficiency(A_t, B_t)
print(f"Cross statistical inefficiency: {g_cross:.2f}")
```

--------------------------------

### Compute Free Energy Differences in pymbar

Source: https://github.com/choderalab/pymbar/blob/main/docs/moving_from_pymbar3.rst

Demonstrates the change in how free energy differences are computed and accessed in pymbar. Version 3 returned a tuple, while version 4 returns a dictionary with keys like 'Delta_f' and 'dDelta_f'. This change standardizes the output format across estimators.

```python
mbar = MBAR(u_kn, N_k)
results, errors = mbar.getFreeEnergyDifferences()
print(results[0])
print(errors[0])
```

```python
mbar = MBAR(u_kn, N_k)
results = mbar.compute_free_energy_differences()
print(results['Delta_f'])
print(results['dDelta_f'])
```

--------------------------------

### Subsample Correlated Data Directly (Python)

Source: https://github.com/choderalab/pymbar/blob/main/docs/timeseries.rst

Directly subsamples a timeseries to obtain an effectively uncorrelated subset of data, without discarding an initial transient. The statistical inefficiency is computed automatically if not provided.

```python
from pymbar import timeseries
indices = timeseries.subsample_correlated_data(A_t_equil)
A_n = A_t_equil[indices]
```

=== COMPLETE CONTENT === This response contains all available snippets from this library. No additional content exists. Do not make further requests.