🌿 SOC-Ensemble-MATLAB

Soil Organic Carbon Estimation from Hyperspectral Reflectance

using Multi-Model Ensemble Machine Learning

Can machine learning read the carbon story hidden in light? This project answers yes — by combining six ML models, multi-transform spectral preprocessing, and auxiliary soil parameters to estimate Soil Organic Carbon from Vis-NIR hyperspectral data with R² = 0.785 and RPD = 2.16.

📊 Results · 🔬 Methodology · 🚀 Quick Start · 📂 Structure · ⚙️ Requirements

---

🎯 Why This Project?

Soil Organic Carbon (SOC) is one of the most critical indicators of soil health, agricultural productivity, and climate change mitigation. Traditional lab-based SOC measurement is:

• 🐢 Slow — weeks of laboratory processing
• 💸 Expensive — wet chemistry reagents, specialized equipment
• 📍 Point-based — cannot scale to field or regional level

Hyperspectral remote sensing changes the game. A single scan across 400–2400 nm captures the molecular fingerprints of carbon bonds, iron oxides, clay minerals, and organic matter — in seconds.

This pipeline builds a production-ready ensemble model that learns from spectral signatures to predict SOC with accuracy comparable to lab methods, at a fraction of the cost and time.

---

📊 Results

🏆 Model Performance Summary

Rank	Model	R²	RMSE (g/kg)	NRMSE (%)	RPD	NS
🥇	Weighted Ensemble	0.7848	0.5975	16.28	2.163	0.78
🥈	Stacking Ensemble	0.7602	0.6307	17.44	2.049	0.75
🥉	Gaussian Process Regression	0.7546	0.6380	18.08	2.025	0.72
4	PLSR (12 components)	0.7202	0.6813	19.67	1.897	0.68
5	Gradient Boosting	0.7134	0.6895	19.43	1.874	0.69
6	Random Forest (tuned)	0.7068	0.6973	18.92	1.853	0.70

RPD interpretation:

• RPD < 1.5 → Poor predictive ability
• 1.5 ≤ RPD < 2.0 → Moderate, useful for screening
• RPD ≥ 2.0 → Good, suitable for quantitative prediction ✅

📈 Visual Results

All 6 Models — Calibration vs. Validation Scatter Plots

PLSR Detail — 12 Latent Components

Spectral Transfer Function (STF) — Spatial Generalization Test

Mean Spectral Reflectance Across SOC Classes (400–2400 nm)

Notice how higher SOC consistently lowers reflectance across the spectrum — this is the physical basis the models exploit.

---

🔬 Methodology

Pipeline Overview

Raw Spectra (400–2400 nm)
         │
         ▼
┌─────────────────────────────┐
│   Multi-Transform Preprocessing  │
│  T1: SG → 1st deriv → SNV   │
│  T2: SG → 2nd deriv → SNV   │
│  T3: MSC → SG → 1st deriv   │
└────────────┬────────────────┘
             │
             ▼
┌─────────────────────────────┐
│   Feature Fusion             │
│  Spectral + [BD, pH, EC,     │
│   CCE, Clay, Sand, Silt]     │
└────────────┬────────────────┘
             │
             ▼
┌─────────────────────────────┐
│   Recursive Feature Elim.    │
│  RF OOB Importance → Top 15% │
└────────────┬────────────────┘
             │
             ▼
┌─────────────────────────────┐
│   10-Fold Stratified CV      │
│  Sorted by SOC → Balanced    │
└────────────┬────────────────┘
             │
      ┌──────┴───────┐
      ▼              ▼
 Base Models     PLSR Model
 RF / GPR / GB   (Top 20% features,
                  1–30 components CV)
      │
      ▼
┌─────────────────────────────┐
│   Ensemble Layer             │
│  Stacking (GPR meta-learner) │
│  Weighted (fminsearch optim) │
└────────────┬────────────────┘
             │
             ▼
     Final SOC Prediction

1️⃣ Spectral Preprocessing — Triple Transform

The raw reflectance is transformed three ways to capture different aspects of the SOC signal:

Transform	Steps	Purpose
T1	SG filter (order=2, win=11) → 1st derivative → SNV	Removes baseline shifts, highlights sharp absorption features
T2	SG filter → 2nd derivative → SNV	Resolves overlapping peaks, enhances subtle bands
T3	MSC → SG filter → 1st derivative	Corrects particle size scattering effects

All three are concatenated into a single feature matrix, tripling the spectral information available to downstream models.

2️⃣ Feature Fusion — Spectral + Soil Physics

Soil carbon doesn't exist in isolation. Bulk density, texture, pH, and carbonate content all modulate the spectral signal. Seven auxiliary parameters are Z-score normalized and appended to the spectral feature matrix:

Feature Matrix = [Spectral Transforms | BD | pH | EC | CCE | Clay | Sand | Silt]

This fusion consistently improves R² by ~0.03–0.05 across all models compared to spectral-only input.

3️⃣ Feature Selection — Recursive Feature Elimination

A preliminary Random Forest (100 trees) computes OOB permutation importance — how much each feature degrades prediction accuracy when randomly shuffled. Features are ranked and the top fraction retained:

• Tree-based models (RF, GPR, GB, Stacking, Weighted): Top 15% (≤200 features)
• PLSR: Top 20% (≤300 features) — preserves more collinear spectral variance

4️⃣ Cross-Validation — Stratified 10-Fold

Standard random k-fold creates folds with unequal SOC distributions, making high-SOC samples systematically under-represented in some folds. This pipeline uses stratified assignment:

1. Sort all samples by SOC value
2. Assign fold labels cyclically: sample 1→fold 1, sample 2→fold 2, ..., sample 11→fold 1, ...
3. Every fold sees the full SOC range from 0.1 to 7.7 g/kg

This eliminates fold-variance artifacts and produces more reliable generalization estimates.

5️⃣ Hyperparameter Tuning — Grid Search for RF

Random Forest is tuned via exhaustive grid search with 5-fold inner CV:

n_trees  ∈ {50, 100, 150, 200}
min_leaf ∈ {5, 8, 10, 12, 15}
→ 20 configurations tested
→ Best config selected by minimum CV-RMSE

6️⃣ Ensemble Strategies

Stacking:

• Base models (RF, GPR, GB) generate out-of-fold predictions
• A GPR meta-learner is trained on these predictions as new features
• Final prediction is the meta-learner's output

Weighted Ensemble:

• Combine all 5 model predictions with learned weights
• Weights optimized via fminsearch minimizing validation RMSE
• Constrained: all weights ≥ 0, sum = 1

---

🔑 Key Spectral Bands

RFE importance analysis consistently highlighted these wavelength regions across all models:

Wavelength (nm)	Physical Assignment	Link to SOC
~495	Iron oxide absorption edge	Indirect — Fe-OM complexes
~730	Red-edge transition	Organic matter chromophores
~970	O-H stretch 2nd overtone	Soil moisture (SOC correlate)
~1416	C-H stretch 1st overtone	Direct organic matter signature
~2200–2400	Al-OH, clay mineral bands	Clay-SOC co-variance

---

📋 Dataset Description

Parameter	Value
Total samples	154
Spectral range	400–2400 nm
Spectral resolution	~1 nm (2001 bands)
SOC range	0.1 – 7.7 g/kg
Low SOC class	0.1–0.7 g/kg (n = 54)
Medium SOC class	0.8–1.4 g/kg (n = 51)
High SOC class	1.6–7.7 g/kg (n = 49)
Auxiliary parameters	7 (BD, pH, EC, CCE, Clay, Sand, Silt)
Train/Test strategy	10-fold stratified CV (no fixed hold-out)

---

⚙️ Requirements

MATLAB Toolboxes

Toolbox	Used For
Statistics and Machine Learning Toolbox	TreeBagger, fitrgp, fitrensemble, plsregress
Signal Processing Toolbox	sgolayfilt
MATLAB core	All matrix operations, fminsearch, polyfit

Tested on MATLAB R2022a and R2023b. Earlier versions may lack some fitrgp options.

Data Files

File	Fields Required
matched_data.mat	wavelengths (1×2001), matched_spectra (2001×n), matched_soc (n×1)
final_data.xlsx	Columns 5–11: BD, pH, EC, CCE, Clay, Sand, Silt

⚠️ Raw field data is not included due to data ownership agreements. To use your own data, prepare the .mat and .xlsx files in the format described above.

---

🚀 Quick Start

%% Step 1 — Place data files in the project folder
%  matched_data.mat
%  final_data.xlsx

%% Step 2 — Open the Live Script
open('soc_ensemble_main.mlx')

%% Step 3 — Update paths at the top of the script
data      = load('matched_data.mat');        % ← update if needed
file_path = 'final_data.xlsx';               % ← update if needed
sheet_name = 'Sheet1';                       % ← match your sheet name

%% Step 4 — Run All (Ctrl + Enter)
%  Expected runtime: ~15–30 min depending on hardware
%  (GPR and grid search are the bottlenecks)

Expected Console Output

=== Loading Data ===
Loaded: 154 samples × 2001 wavelengths
=== Loading Auxiliary Soil Parameters ===
✓ Auxiliary parameters loaded: 154 samples
=== Wavelength Selection ===
Wavelengths: 2001 → 2001
=== Multi-Transform Preprocessing ===
=== Recursive Feature Elimination ===
Features selected: 200
=== Model 1: Grid Search RF ===
  Testing: Trees=50, MinLeaf=5... RMSE=0.7214
  ...
Best RF: Trees=150, MinLeaf=8, R²=0.7068, RMSE=0.6973
=== Model 2: Gaussian Process Regression ===
...
=== FINAL ENSEMBLE RESULTS ===
Weighted Ensemble: R²=0.7848, RMSE=0.5975, RPD=2.163

---

📂 Repository Structure

SOC-Ensemble-MATLAB/
│
├── 📄 soc_ensemble_main.mlx              # Main MATLAB Live Script
│
├── 🖼️ Screenshot 2026-07-12 202023.png   # All 6 models — calibration & validation
├── 🖼️ Screenshot 2026-07-12 202038.png   # PLSR scatter (12 components)
├── 🖼️ Screenshot 2026-07-12 202052.png   # STF spatial validation
├── 🖼️ Screenshot 2026-07-12 202105.png   # Mean reflectance by SOC class
│
├── 📄 README.md
└── 📄 LICENSE

---

🧠 Technical Highlights

• ✅ Multi-transform preprocessing — three spectral transforms concatenated for richer features
• ✅ Feature fusion — spectral + 7 auxiliary soil parameters
• ✅ Two-stage RFE — different feature budgets for linear vs. nonlinear models
• ✅ Stratified CV — SOC-balanced folds for reliable generalization estimates
• ✅ Grid search hyperparameter tuning — 20 RF configurations evaluated
• ✅ Two ensemble strategies — stacking and optimized weighting
• ✅ Full metrics suite — R², RMSE, NRMSE, RPD, NS reported for all models

---

📄 License

This project is licensed under the MIT License — see LICENSE for full terms. You are free to use, modify, and distribute this code with attribution.

---

👤 Author

Sheyda Asadi

🔗 github.com/sheyda2021

---

If this project was useful or interesting to you, a ⭐ star on GitHub goes a long way — it helps others find the work and means a lot!