GEC Monopore Characteristic Function

8.2. GEC Monopore Characteristic Function#

This notebook derives the GEC (Giddings–Eyring–Carmichael) monopore characteristic function from first principles using SymPy.

Goal: Show how the CF \(\phi_{t_S}(\omega) = \exp\left[\bar{n}\left(\frac{1}{1-i\omega\bar{\tau}} - 1\right)\right]\) emerges from:

  1. Exponential egress time distribution

  2. Poisson number of adsorptions

  3. Independence assumption \(r_M \perp \{\tau_{S,i}\}\)

Attribution

This educational notebook was developed with assistance from GitHub Copilot (Claude Sonnet 4.5), December 2024. The step-by-step derivation and symbolic verification code were collaboratively created to demonstrate the mathematical foundations of the GEC model for chromatography applications.

8.2.1. Step 1: CF of Exponential Distribution#

For a single egress time \(\tau \sim \text{Exp}(1/\bar{\tau})\):

\[\phi_\tau(\omega) = \mathbb{E}[e^{i\omega\tau}] = \int_0^\infty e^{i\omega t} \frac{1}{\bar{\tau}} e^{-t/\bar{\tau}} dt = \frac{1}{1 - i\omega\bar{\tau}}\]
from sympy import symbols, exp, factorial, simplify, I, oo, Sum
import sympy as sp

# Define symbols
w, n, t_bar, n_bar = symbols('w n t_bar n_bar', real=True, positive=True)

8.2.2. Step 2: CF of Sum of n i.i.d. Exponentials#

For \(S_n = \sum_{i=1}^n \tau_i\) where \(\tau_i\) are i.i.d. exponential:

\[\phi_{S_n}(\omega) = [\phi_\tau(\omega)]^n = \left(\frac{1}{1-i\omega\bar{\tau}}\right)^n\]

This uses independence: \(\mathbb{E}[e^{i\omega(X+Y)}] = \mathbb{E}[e^{i\omega X}]\mathbb{E}[e^{i\omega Y}]\)

# CF of single exponential random variable
phi_tau = 1 / (1 - I*w*t_bar)
print("CF of single egress time:")
phi_tau

CF of single egress time:
\[\displaystyle \frac{1}{- i t_{bar} w + 1}\]

8.2.3. Step 3: Poisson Distribution#

For \(r_M \sim \text{Poisson}(\bar{n})\):

\[P(r_M = n) = \frac{e^{-\bar{n}}\bar{n}^n}{n!}\]
# CF of sum of n i.i.d. exponentials
phi_sum_n = phi_tau**n
print("CF of sum of n egress times:")
phi_sum_n

CF of sum of n egress times:
\[\displaystyle \left(- i t_{bar} w + 1\right)^{- n}\]

8.2.4. Step 4: Law of Total Expectation (Requires Independence!)#

For random sum \(t_S = \sum_{i=1}^{r_M} \tau_i\) where \(r_M \perp \{\tau_i\}\):

\[\phi_{t_S}(\omega) = \mathbb{E}_{r_M}\left[\mathbb{E}[e^{i\omega t_S} | r_M]\right] = \sum_{n=0}^\infty P(r_M=n) \cdot [\phi_\tau(\omega)]^n\]

Critical: This step requires independence between \(r_M\) and \(\{\tau_i\}\)!

# Poisson probability mass function
P_n = exp(-n_bar) * n_bar**n / factorial(n)
print("Poisson PMF P(r_M = n):")
P_n

Poisson PMF P(r_M = n):
\[\displaystyle \frac{n_{bar}^{n} e^{- n_{bar}}}{n!}\]

8.2.5. Step 5: Evaluate the Sum (Poisson Generating Function)#

We need to evaluate:

\[\sum_{n=0}^\infty \frac{e^{-\bar{n}}\bar{n}^n}{n!} \cdot z^n \quad \text{where } z = \phi_\tau(\omega)\]

This is the probability generating function of Poisson:

\[\mathbb{E}[z^{r_M}] = e^{\bar{n}(z-1)}\]
# Compound Poisson CF: sum over all possible n values
# φ_tS(ω) = Σ P(n) * [φ_τ(ω)]^n
summand = P_n * phi_sum_n
print("Summand P(r_M=n) * [φ_τ(ω)]^n:")
summand

Summand P(r_M=n) * [φ_τ(ω)]^n:
\[\displaystyle \frac{n_{bar}^{n} \left(- i t_{bar} w + 1\right)^{- n} e^{- n_{bar}}}{n!}\]

8.2.6. Step 6: Final Result (Dondi Eq. 43)#

Substituting \(z = \frac{1}{1-i\omega\bar{\tau}}\) into \(e^{\bar{n}(z-1)}\):

\[\phi_{t_S}(\omega) = \exp\left[\bar{n}\left(\frac{1}{1-i\omega\bar{\tau}} - 1\right)\right]\]

This can also be written as:

\[\phi_{t_S}(\omega) = \exp\left[\frac{\bar{n} \cdot i\omega\bar{\tau}}{1-i\omega\bar{\tau}}\right]\]
# The sum Σ (e^(-n_bar) * n_bar^n / n!) * z^n = exp(n_bar*(z-1))
# We substitute z = phi_tau
z = symbols('z')

# Poisson generating function
pgf = exp(n_bar * (z - 1))
print("Poisson generating function E[z^r_M]:")
display(pgf)

# Substitute z = phi_tau to get GEC CF
gec_cf_derived = pgf.subs(z, phi_tau)
print("\nGEC CF (substitute z = φ_τ(ω)):")
gec_cf_derived

Poisson generating function E[z^r_M]:
\[\displaystyle e^{n_{bar} \left(z - 1\right)}\]
GEC CF (substitute z = φ_τ(ω)):
\[\displaystyle e^{n_{bar} \left(-1 + \frac{1}{- i t_{bar} w + 1}\right)}\]

8.2.7. Verification: Check Against Original Form#

Now we verify this matches the form used in later cells.

# Simplify the exponent
exponent = simplify(phi_tau - 1)
print("Exponent φ_τ(ω) - 1:")
display(exponent)

# Full GEC CF
gec_cf_final = simplify(gec_cf_derived)
print("\nFinal GEC monopore CF:")
gec_cf_final

Exponent φ_τ(ω) - 1:
\[\displaystyle - \frac{t_{bar} w}{t_{bar} w + i}\]
Final GEC monopore CF:
\[\displaystyle e^{- \frac{i n_{bar} t_{bar} w}{i t_{bar} w - 1}}\]
from sympy import symbols, exp, simplify, I, cancel
w = symbols('w')
n1, t1 = symbols('n1, t1', positive=True, real=True)

# Original form from Dondi Eq. 43
gec_monopore_cf_original = exp(n1*(1/(1 - I*w*t1) - 1))

print("Original form (Dondi Eq. 43):")
display(gec_monopore_cf_original)

# Manually simplify the exponent to show equivalence
# 1/(1 - iωτ) - 1 = [1 - (1 - iωτ)]/(1 - iωτ) = iωτ/(1 - iωτ)
exponent_original = n1*(1/(1 - I*w*t1) - 1)
exponent_simplified = simplify(exponent_original)

print("\nOriginal exponent:")
display(exponent_original)
print("\nSimplified exponent:")
display(exponent_simplified)

# This should be: n1*I*w*t1/(1 - I*w*t1)
# Verify by manual calculation
numerator = 1 - (1 - I*w*t1)
denominator = 1 - I*w*t1
manual_simplified = n1 * numerator / denominator
manual_simplified = simplify(manual_simplified)

print("\nManual simplification [1 - (1 - iωτ)]/(1 - iωτ):")
display(manual_simplified)

# Alternative form
gec_alternative = exp(n1*I*w*t1/(1 - I*w*t1))
print("\nAlternative form exp[n̄·iω·τ̄/(1 - iω·τ̄)]:")
display(gec_alternative)

# Verify they're the same by checking the difference
print("\nDifference between original and alternative:")
display(simplify(gec_monopore_cf_original - gec_alternative))

print("\n✓ Both forms are mathematically identical!")
print("  Form 1: exp[n̄(1/(1-iωτ̄) - 1)]")
print("  Form 2: exp[n̄·iωτ̄/(1-iωτ̄)]")

Original form (Dondi Eq. 43):
\[\displaystyle e^{n_{1} \left(-1 + \frac{1}{- i t_{1} w + 1}\right)}\]
Original exponent:
\[\displaystyle n_{1} \left(-1 + \frac{1}{- i t_{1} w + 1}\right)\]
Simplified exponent:
\[\displaystyle - \frac{n_{1} t_{1} w}{t_{1} w + i}\]
Manual simplification [1 - (1 - iωτ)]/(1 - iωτ):
\[\displaystyle - \frac{n_{1} t_{1} w}{t_{1} w + i}\]
Alternative form exp[n̄·iω·τ̄/(1 - iω·τ̄)]:
\[\displaystyle e^{\frac{i n_{1} t_{1} w}{- i t_{1} w + 1}}\]
Difference between original and alternative:
\[\displaystyle 0\]
✓ Both forms are mathematically identical!
  Form 1: exp[n̄(1/(1-iωτ̄) - 1)]
  Form 2: exp[n̄·iωτ̄/(1-iωτ̄)]

8.2.8. Moments from Characteristic Functions#

The characteristic function \(\phi(\omega)\) encodes all statistical properties of a distribution. We can extract moments using derivatives:

Raw Moments (about the origin)#

The \(k\)-th raw moment is \(\mu'_k = \mathbb{E}[X^k]\), obtained via:

\[\mu'_k = (-i)^k \left.\frac{d^k \phi(\omega)}{d\omega^k}\right|_{\omega=0}\]

  • \(k=1\): Mean (first raw moment)

  • \(k=2\): Second raw moment \(\mathbb{E}[X^2]\)

Central Moments (about the mean)#

The \(k\)-th central moment is \(\mu_k = \mathbb{E}[(X-\mu)^k]\) where \(\mu = \mathbb{E}[X]\):

  • \(k=2\): Variance \(\sigma^2\)

  • \(k=3\): Skewness component (normalized: \(\gamma_1 = \mu_3/\sigma^3\))

  • \(k=4\): Kurtosis component (normalized: \(\kappa = \mu_4/\sigma^4\))

Method: Central moments can be computed from the shifted CF:

\[\tilde{\phi}(\omega) = \phi(\omega) \cdot e^{-i\omega\mu} \quad \Rightarrow \quad \mu_k = (-i)^k \left.\frac{d^k \tilde{\phi}(\omega)}{d\omega^k}\right|_{\omega=0}\]

This is implemented below by dividing the CF by \(e^{i\omega\mu}\) before differentiation.

8.2.9. How to use SymPy to symbolically compute moments#

We can extract statistical moments from the characteristic function using SymPy as follows.

from sympy import symbols, exp, diff, simplify, I
w = symbols('w')

def raw_moment(cf, k):
    return simplify((-I)**k * diff(cf, w, k).subs(dict(w=0)))

def central_moment(cf, k):
    m = raw_moment(cf, 1)
    return simplify(raw_moment(cf/exp(I*w*m),k))

n1, t1 = symbols('n1, t1')
gec_monopore_cf = exp(n1*(1/(1 - I*w*t1) - 1))
gec_monopore_cf
\[\displaystyle e^{n_{1} \left(-1 + \frac{1}{- i t_{1} w + 1}\right)}\]
raw_moment(gec_monopore_cf, 1)
\[\displaystyle n_{1} t_{1}\]
central_moment(gec_monopore_cf, 2)
\[\displaystyle 2 n_{1} t_{1}^{2}\]
central_moment(gec_monopore_cf, 3)
\[\displaystyle 6 n_{1} t_{1}^{3}\]
Contents