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}")