Skip to content

Chemical Composition Analysis

Chemical composition analysis workflow, including mass reconstruction, ion balance, hygroscopicity calculation, and gas-particle partitioning ratios.

Data Preparation

from datetime import datetime
from pathlib import Path
import pandas as pd
from AeroViz import (
    RawDataReader,
    reconstruct_mass,
    volume_ri,
    kappa,
    growth_factor,
    partition_ratios,
    split_oc_ec,
)

# Read ion data
igac = RawDataReader(
    instrument='IGAC',
    path=Path('./data/igac'),
    start=datetime(2024, 1, 1),
    end=datetime(2024, 3, 31)
)

# Read carbon component data
ocec = RawDataReader(
    instrument='OCEC',
    path=Path('./data/ocec'),
    start=datetime(2024, 1, 1),
    end=datetime(2024, 3, 31)
)

# Read elemental data (Xact 625i)
xact = RawDataReader(
    instrument='Xact',
    path=Path('./data/xrf'),
    start=datetime(2024, 1, 1),
    end=datetime(2024, 3, 31)
)

# Merge data
df_chem = pd.concat([igac, ocec, xact], axis=1)

Mass Reconstruction

Basic Reconstruction

# Perform mass reconstruction
result = reconstruct_mass(df_chem, df_ref=df_chem[['PM25']])

# Main component masses
df_mass = result['mass']
print(df_mass.columns)
# ['AS', 'AN', 'OM', 'Soil', 'SS', 'EC', 'total']
# 'total' = reconstructed PM mass (sum of species) — there is no 'PM25_rc' column

# Ammonium status: a DataFrame with 'ratio' and 'status' columns
nh4_status = result['NH4_status']
print(nh4_status['status'].value_counts())
# Enough        75
# Deficiency    15

Closure Check

# Calculate closure ('total' is the reconstructed mass)
closure = df_mass['total'] / df_chem['PM25'] * 100

print(f"Mean closure: {closure.mean():.1f}%")
print(f"Std closure: {closure.std():.1f}%")

# Closure distribution
import matplotlib.pyplot as plt
plt.hist(closure, bins=20, edgecolor='black')
plt.xlabel('Closure (%)')
plt.ylabel('Frequency')
plt.axvline(100, color='r', linestyle='--', label='100%')
plt.legend()
plt.show()

Component Contributions

# Calculate contribution ratios for each component
components = ['AS', 'AN', 'OM', 'EC', 'Soil', 'SS']
contributions = df_mass[components].div(df_mass['total'], axis=0) * 100

print("Average contribution (%):")
print(contributions.mean())
# AS      25.3
# AN      18.2
# OM      35.1
# EC       8.5
# Soil     7.2
# SS       5.7

Volume and Refractive Index Calculation

# Reconstruction already produced volumes; feed them to volume_ri
df_volume = result['volume']
print(df_volume.columns)
# ['AS_volume', 'AN_volume', 'OM_volume', 'EC_volume', 'Soil_volume', 'SS_volume', 'total_dry']

# Volume-weighted refractive index (dry + ambient if ALWC provided)
df_RI = volume_ri(df_volume, df_alwc=df_alwc)
print(f"Mean n_dry: {df_RI['n_dry'].mean():.3f}")
print(f"Mean k_dry: {df_RI['k_dry'].mean():.4f}")

Hygroscopicity (kappa) Calculation

import pandas as pd

# Growth factor from volume + ALWC
df_gRH = growth_factor(df_volume, df_alwc)        # column: gRH
print(f"Mean gRH: {df_gRH['gRH'].mean():.2f}")

# kappa: needs gRH, ambient temperature (AT, °C), and RH (%)
df_kappa_input = pd.concat([df_gRH, met_data[['AT', 'RH']]], axis=1)
df_kappa = kappa(df_kappa_input, diameter=0.5)    # column: kappa_chem
print(f"Mean kappa: {df_kappa['kappa_chem'].mean():.3f}")

kappa vs Composition Relationship

# Analyze relationship between kappa and chemical composition
fig, axes = plt.subplots(1, 3, figsize=(12, 4))

# kappa vs SIA ratio
sia_ratio = (df_mass['AS'] + df_mass['AN']) / df_mass['total']
axes[0].scatter(sia_ratio, df_kappa['kappa_chem'], alpha=0.5)
axes[0].set_xlabel('SIA / PM2.5')
axes[0].set_ylabel('kappa')

# kappa vs OM ratio
om_ratio = df_mass['OM'] / df_mass['total']
axes[1].scatter(om_ratio, df_kappa['kappa_chem'], alpha=0.5)
axes[1].set_xlabel('OM / PM2.5')
axes[1].set_ylabel('kappa')

# kappa vs RH
axes[2].scatter(met_data['RH'], df_kappa['kappa_chem'], alpha=0.5)
axes[2].set_xlabel('RH (%)')
axes[2].set_ylabel('kappa')

plt.tight_layout()
plt.show()

Gas-Particle Partitioning Ratios

# Requires gas data + temperature
df_gas = gas_data[['SO2', 'NO2', 'HNO3', 'NH3']]
df_combined = pd.concat([df_chem, df_gas, met_data[['temp']]], axis=1)

# Calculate partitioning ratios
partition = partition_ratios(df_combined)

print(partition.columns)
# ['SOR', 'NOR', 'NOR_2', 'NTR', 'epls_SO42-', 'epls_NO3-', 'epls_NH4+', 'epls_Cl-']

Partitioning Ratio Description

Indicator Formula Meaning
SOR SO4^2-/(SO4^2-+SO2) Sulfate conversion degree
NOR NO3-/(NO3-+NO2) Nitrate conversion degree
NTR NO3-/(NO3-+HNO3) Nitrate partitioning
epsilon_ite NO3-/(NO3-+Cl-) Nitrate vs chloride
# Analyze diurnal variation of SOR and NOR
hourly_sor = partition['SOR'].groupby(partition.index.hour).mean()
hourly_nor = partition['NOR'].groupby(partition.index.hour).mean()

fig, ax = plt.subplots()
ax.plot(hourly_sor.index, hourly_sor.values, 'b-o', label='SOR')
ax.plot(hourly_nor.index, hourly_nor.values, 'r-o', label='NOR')
ax.set_xlabel('Hour')
ax.set_ylabel('Ratio')
ax.legend()
plt.show()

OC/EC Split (POC / SOC)

# Requires the OC/EC analyzer level-result columns
# (OC1-OC4, PC, Thermal_OC/EC, Optical_OC/EC, Sample_Volume).
ocec_result = split_oc_ec(df_lcres, df_mass=df_chem[['PM25']])

# Concatenated OC/EC + POC/SOC/WSOC/WISOC + status flags
basic = ocec_result['basic']
print(f"Mean OC/EC: {basic['OC/EC'].mean():.2f}")
print(f"Mean SOC: {basic['SOC'].mean():.2f} ug/m3")

# Per-species PM ratios
ratios = ocec_result['ratio']

Complete Analysis Script

from datetime import datetime
from pathlib import Path
import pandas as pd
import matplotlib.pyplot as plt
from AeroViz import (
    RawDataReader,
    reconstruct_mass,
    volume_ri,
    kappa,
    growth_factor,
)

# 1. Read data (pass dates as start=/end= keywords)
igac = RawDataReader('IGAC', Path('./data/igac'),
                     start='2024-01-01', end='2024-03-31', mean_freq='1h')
ocec = RawDataReader('OCEC', Path('./data/ocec'),
                     start='2024-01-01', end='2024-03-31', mean_freq='1h')
df_chem = pd.concat([igac, ocec], axis=1)
met    = pd.read_csv('met.csv', index_col='time', parse_dates=True)
df_alwc = met[['ALWC']]                                  # from ISOROPIA, e.g.

# 2. Mass reconstruction
mass_result = reconstruct_mass(df_chem, df_ref=df_chem[['PM25']])
df_mass   = mass_result['mass']
df_volume = mass_result['volume']

# 3. Volume-weighted refractive index
df_RI = volume_ri(df_volume, df_alwc=df_alwc)

# 4. Hygroscopicity
df_gRH   = growth_factor(df_volume, df_alwc)
df_kappa = kappa(
    pd.concat([df_gRH, met[['AT', 'RH']]], axis=1),
    diameter=0.5,
)

# 5. Output results summary
print("=== Chemical Analysis Summary ===")
print(f"PM2.5: {df_chem['PM25'].mean():.1f} +/- {df_chem['PM25'].std():.1f} ug/m3")
print(f"Closure: {(df_mass['total']/df_chem['PM25']*100).mean():.1f}%")
print(f"RI(dry): {df_RI['n_dry'].mean():.3f} + {df_RI['k_dry'].mean():.4f}i")
print(f"kappa: {df_kappa['kappa_chem'].mean():.3f}")