Hyperspectral Band Selection for Ground Fuel Classification for Prescribed Fires

1. Introduction

Hyperspectral imaging (HSI) has emerged as a powerful remote sensing tool, enabling fine-grained discrimination of materials, vegetation health assessment, and monitoring of environmental disturbances due to its hundreds of contiguous narrow spectral bands [1,2,3]. The rich spectral information supports applications ranging from land cover classification [4] and crop mapping [5] to soil analysis [6] and quality inspection in agriculture [7]. However, the high dimensionality of HSI data introduces significant challenges, including the “curse of dimensionality,” spectral redundancy, computational complexity, and reduced model generalization in machine learning pipelines [8,9].

Band selection techniques address these issues by identifying a compact, informative subset of bands that preserve discriminative power while minimizing redundancy and noise [1,10]. Traditional approaches, such as Principal Component Analysis (PCA), rely on variance maximization but often fail to capture task-specific discriminability [3]. More advanced methods incorporate spatial-spectral information [11], attention mechanisms [12], or graph-based structures [13]. Recently, deep reinforcement learning (DRL) has shown promise for adaptive, task-driven band selection in HSI classification, outperforming conventional techniques in semisupervised and unsupervised scenarios [14,15,16,17].

These challenges are particularly acute in wildfire science and prescribed fire management, where spectral responses are influenced by burn severity, char deposition, canopy structure, and early vegetation recovery [18,19]. Wildfires pose a growing threat to ecosystems, communities, and economies in the western United States, with Montana experiencing particularly pronounced impacts due to its vast forested landscapes, dry continental climate, and increasing human development in the wildland-urban interface. Over the past two decades, Montana has seen a dramatic rise in wildfire activity: average annual acres burned have exceeded 300,000 since 2000–more than ten times the pre-1950 average of less than 30,000 acres–driven by warmer temperatures, prolonged droughts, fuel accumulation from historical fire suppression, and climate change-induced aridity [20,21]. Recent seasons illustrate variability but underscore the trend: in 2024, over 352,000 acres burned, while 2025 saw a milder year with approximately 76,000 acres (one of the fourth-lowest totals in the last 15 years), attributed to favorable weather but still highlighting ongoing risks [22,23].

Prescribed fires–carefully planned, low-intensity burns–serve as a critical proactive tool to mitigate these threats by reducing hazardous fuel loads, restoring ecological processes, improving wildlife habitat, recycling nutrients, and promoting resilient forests. In Montana, where historical natural fire regimes burned tens to hundreds of thousands of acres annually before suppression, prescribed burning helps counteract fuel buildup and supports forest health. State and federal programs (e.g., DNRC, USFS Helena-Lewis and Clark National Forest) have expanded efforts, with recent accomplishments including tens of thousands of acres treated annually and proposals to scale to approximately 40,000 acres per year through 2045 on national forests alone [24,25,26]. Despite these advances, challenges persist, including limited burn windows due to air quality regulations, smoke concerns, and the need for precise monitoring of fuel types and post-burn recovery.

Accurate ground fuel classification–distinguishing live vegetation, dead biomass, and species-specific conditions–is essential for planning effective prescribed burns, predicting fire behavior, and evaluating ecological outcomes. Hyperspectral remote sensing, with its fine spectral resolution, offers significant potential for this task, yet real-world applications in fire-affected Montana forests remain underexplored due to data scarcity and dimensionality issues [27,28,29].

Despite advances in HSI band selection for specific tasks [30,31], few studies integrate multiple strategies (e.g., statistical, edge-aware, attention-guided, and reinforcement learning-based) in a unified pipeline, especially for prescribed fire applications. Prior work has explored attention and edge-aware methods for burned vegetation classification [32] and self-supervised/DRL approaches for prescribed burn impact analysis [33], highlighting the potential of learning-based techniques in fire-affected scenes. However, systematic comparisons across benchmarks and real-world Unmanned Aircraft System (UAS) datasets remain limited, particularly for ground fuel mapping in controlled burns.

This study addresses these gaps by presenting a comprehensive evaluation of four band-selection strategies–Principal Component Analysis (PCA), Spatial-Spectral Edge Preservation (SSEP), Spectral-Redundancy Penalized Attention (SRPA), and Deep Reinforcement Learning (DRL)–combined with classical (Random Forest, Support Vector Machines, K-Nearest Neighbors) and deep (3D Convolutional Neural Networks) classifiers. We apply a modular, reproducible pipeline to benchmark datasets (Indian Pines, Pavia University, Salinas, Botswana, Kennedy Space Center) and a novel VNIR hyperspectral dataset collected via UAS over prescribed burn sites at the Lubrecht Experimental Forest, Montana, USA. This dataset captures pre-burn conditions in a controlled thin-burn plot, providing a unique real-world testbed for understanding band-selection behavior in fire-prone forest environments.

By emphasizing exploratory data analysis (EDA), preprocessing, patch-based modeling, and cross-dataset validation, our work identifies compact spectral subsets that enhance classification accuracy while reducing computational demands. The inclusion of the Montana dataset bridges the gap between controlled benchmarks and operational remote sensing in wildfire contexts, contributing a generalizable framework for HSI dimensionality reduction and classification with applications in prescribed fire monitoring, vegetation recovery assessment, and fuel load management.

2. Materials and Methods

The goal of this study is to systematically evaluate a unified and reproducible hyperspectral image (HSI) classification pipeline across both benchmark airborne datasets and a complex UAS-based visible-near-infrared (VNIR) dataset captured over prescribed burn sites in the Lubrecht Experimental Forest, a large outdoor forest research laboratory located in the Blackfoot River drainage, about 30 miles northeast of Missoula, Montana. To accomplish this, we designed a modular workflow consisting of four tightly integrated stages illustrated in process diagram Figure 1: (1) dataset preparation and exploratory data analysis (EDA), (2) noise detection and data cleaning, (3) band-selection using classical and deep-learning-based feature reduction techniques, and (4) classification using both traditional machine learning models and deep learning architectures. All components were implemented in a consistent, dataset-agnostic fashion so that the same pipeline could be applied to Indian Pines, Pavia University, Salinas, Botswana, KSC, and the Montana VNIR dataset without modification.

This section describes the materials, datasets, preprocessing steps, algorithms, model architectures, and evaluation design used in the study, presenting each part as a connected narrative that reflects the logical progression of the project.

2.1. Data

Hyperspectral imagery consists of hundreds of narrow and contiguous spectral bands, providing rich spectral information for material discrimination. However, this high spectral resolution also introduces significant redundancy and the well-known Hughes phenomenon (also known as the "curse of dimensionality"), motivating the need for dimensionality reduction and band selection prior to classification. Consequently, standardized preprocessing and band selection are widely adopted in hyperspectral image analysis pipelines [1,2,3].

We evaluate a unified hyperspectral classification framework using a diverse collection of datasets that span controlled benchmark scenes and a complex real-world unmanned aircraft system (UAS)–based acquisition. Specifically, experiments were conducted on five widely used benchmark hyperspectral datasets—Indian Pines, Pavia University, Salinas Valley, Botswana, and Kennedy Space Center [34]—as well as a custom visible–near-infrared (VNIR) hyperspectral dataset collected over prescribed burn sites in Montana, USA. The benchmark datasets represent a range of land-cover types, spatial resolutions, sensors, and class distributions, and have been extensively used in prior hyperspectral classification and band-selection studies, enabling meaningful comparison with existing literature [8,9,35].

A summary of the benchmark datasets, including sensor type, spatial dimensions, number of spectral bands before and after sensor-documented noisy-band removal, number of classes, and spatial resolution, is provided in Table 1. The cleaned benchmark cubes served as inputs to the subsequent exploratory data analysis. Rather than treating each benchmark dataset independently, all were processed using the same unified and dataset-agnostic pipeline to ensure consistency across preprocessing, band selection, and classification stages. This design allows performance differences to be attributed to algorithmic choices rather than dataset-specific handling.

In contrast to the benchmark scenes, the Montana VNIR dataset represents a real-world post-fire forest environment acquired using a UAS platform, characterized by higher spectral noise, severe class imbalance, sparse ground truth, and complex canopy structure. Detailed acquisition parameters, preprocessing steps, and dataset-specific challenges for the Montana dataset are described separately in Section 2.1.1, as they differ substantially from those of the benchmark datasets and require specialized treatment.

2.1.1. Montana

The hyperspectral images used in this study were acquired on 11 May 2024 over a controlled thin-burn plot (330 m × 300 m) in the Lubrecht Experimental Forest, approximately 42 km east of Missoula, Montana, USA. Data collection was conducted using a Headwall co-aligned visible–near infrared (VNIR) and short-wave infrared (SWIR) imaging system, integrated with a FreeFly ALTA X heavy-duty unmanned aircraft system (UAS). Prior to the flight, the imaging system was calibrated using a certified Spectralon white diffuse reflectance standard to determine the optimal exposure time, frame period, and recommended flight speed. Based on this information, a flight plan was designed to maintain a nearly constant altitude relative to the terrain, thereby eliminating the need for terrain correction in the resulting imagery. Magnetometer calibration was also performed to ensure accurate detection of the Earth’s magnetic field for reliable navigation. The precise location of the imaging system was recorded using a Trimble SPS585 global navigation satellite system (GNSS) unit for Real-Time Kinematic (RTK) surveying. Before executing the flight plan, the Inertial Measurement Unit (IMU) was initialized to track the aircraft’s movement and orientation. The hyperspectral imaging system captured 543 spectral band images — 273 in the VNIR range (400–1000 nm) and 270 in the SWIR range (900–2500 nm). The data files acquired were preprocessed using Spectral View, the proprietary software provided by the manufacturer of the hyperspectral imaging system. The raw data, initially recorded in digital numbers (DNs), was corrected for noise and dark current using imagery collected prior to data acquisition. This process converted the DN values into spectral radiance (mW/(cm²·sr·µm)). The radiance data was then transformed into spectral reflectance using a calibrated reference tarp placed within the imaging scene during processing. Finally, the processed images were exported in GeoTIFF format to enable further analysis using standard image processing software.

2.2. Exploratory Data Analysis

Exploratory data analysis (EDA) was conducted prior to model training to examine radiometric distributions, assess spectral redundancy, and identify noisy or unstable bands. EDA plays a critical role in ensuring data integrity and preventing the propagation of artifacts into downstream band selection and classification model [9,27].

2.3. EDA for Benchmark Hyperspectral Datasets

The objective of the benchmark dataset EDA was to standardize all preprocessing steps across the five classical hyperspectral scenes—Indian Pines, Pavia University, Salinas, Botswana, and Kennedy Space Center (KSC)—so that band–selection algorithms could be evaluated under consistent noise conditions. These scenes differ substantially in spatial dimensions, radiometric behavior, and class distributions, making a unified and statistically grounded EDA pipeline necessary. The goal of this pipeline was to verify data integrity, characterize radiometric distributions, identify low–SNR or water–absorption bands, and construct noise–free normalized cubes to serve as inputs for SSEP, SRPA, and DRL–based band selection as well as for classical and deep classifiers.

Each dataset was loaded from its .mat hyperspectral file and associated ground–truth (gt) .mat label mask using dataset–specific variable keys. Upon loading, the benchmark hyperspectral cubes exhibited the following dimensions: Indian Pines (145 × 145 × 200), PaviaU (610 × 340 × 103), Salinas (512 × 217 × 204), Botswana (1476 × 256 × 145) and KSC (512 × 614 × 176). Basic integrity diagnostics included computation of global intensity minima and maxima, band–wise means, and band-wise standard deviations. Botswana displayed reflectance values ranging from 0 to 45106, with an average band-wise mean of 1749.62 and an average standard deviation of 420.48. Indian Pines ranged from 955 to 9604, with an average band–wise mean 2652.39 and an average standard deviation 336.36. KSC ranged from –27 to 1244, with average band–wise mean 93.42 and an average standard deviation of 72.36. PaviaU showed intensities between 0 and 8000, with an average band-wise mean of 1389.12 and an average standard deviation of 716.43. Salinas ranged from –11 to 9207, with corresponding mean band-wise mean of 1196.40 and an average standard deviation of 389.91. These statistics confirmed proper radiometric scaling and the absence of corruption in the raw reflectance values.

Ground-truth label distributions were analyzed to verify class presence and quantify class imbalance. Botswana contained 15 classes with pixel counts ranging from 95 to 314, against 374,608 background pixels. Indian Pines exhibited 16 classes, including extremely sparse categories (e.g., class 1 with only 46 pixels) and dense categories such as class 11 with 2455 pixels. KSC contained 13 classes with class populations ranging from 105 to 927 and background comprising 309,157 pixels. PaviaU included 9 semantic classes, with some highly populated categories (e.g., class 2 with 18,649 pixels) and background containing 164,624 pixels. Salinas exhibited 16 classes, with major categories (e.g., class 8 with 11,271 pixels) and small categories (e.g., class 13 with 916 pixels). These results confirm substantial class imbalance across all benchmark datasets, which should be considered in downstream learning.

To characterize spectral noise, each hyperspectral cube was reshaped into a matrix of size $N\times B$, where $N=H·W$, enabling band–wise analysis across all spatial pixels. For each spectral channel b, we computed its minimum and maximum intensities $(m_b,M_b)$, mean ${\mu }_b$, standard deviation ${\sigma }_b$, dynamic range $R_b=M_b-m_b$, and zero–fraction

$$f_b^{\left(0\right)}={\displaystyle \frac{1}{N}}\sum _{i=1}^N\mathbb{1}\left[X_{i,b}=0\right]$$

(1)

Bands were flagged as suspicious if they exhibited extremely low variance (i.e., ${\sigma }_b$ lying in the lowest 5% percentile of all variances), unusually small dynamic range (i.e., $R_b$ less than ${10}^{-3}$ of the global range), or high zero–fraction (i.e., $f_b^{\left(0\right)}\ge 0.80$). Although none of the benchmark datasets contained high zero–fraction bands, each dataset presented a small set of low–variance bands that aligned with known water–absorption and sensor–instability regions. This EDA-based noisy-band identification was performed in addition to the initial sensor-documented band exclusion summarized in Table 1. The heuristic procedure identified the following candidate noisy bands: Botswana—[111, 138–144]; Indian Pines—[102–104, 143–145, 195, 197–199]; KSC—[0–8]; PaviaU—[0–5]; Salinas—[106–109, 146–148, 198, 201–203]. These empirically identified indices match previously documented noisy–band lists in the hyperspectral imaging literature.

Based on these technical diagnostics and literature support, we removed the identified noisy bands from each dataset, yielding cleaned spectral dimensionalities of 137 bands (Botswana), 190 bands (Indian Pines), 167 bands (KSC), 97 bands (PaviaU), and 193 bands (Salinas). Spatial dimensions were unchanged, ensuring that only the spectral axis was compressed. This removal eliminated known water–absorption wavelengths and low–SNR regions while preserving all meaningful spectral–spatial structure required for downstream tasks.

After noisy–band removal, each dataset was normalized per spectral channel using min–max scaling:

$$X_{i,b}^{\mathrm{norm}}={\displaystyle \frac{X_{i,b}-m_b}{M_b-m_b+\epsilon }}$$

(2)

where $\epsilon$ prevents division by zero. This transformation places all spectral channels in the $[0,1]$ range, eliminates scale disparities among bands, and stabilizes distance metrics and gradients for machine learning models. The cleaned and normalized hyperspectral cubes and associated summary statistics were retained for downstream analysis.

The results of this benchmark EDA demonstrate that the unified pipeline not only confirms the structural and radiometric consistency of all datasets but also quantitatively identifies and removes noisy spectral regions using both heuristic and literature–validated criteria. By producing noise–controlled hyperspectral cubes with well–characterized statistics, the benchmark EDA ensures that subsequent band–selection and classification experiments are conducted under standardized, reproducible, and scientifically rigorous conditions.

2.4. EDA for the Montana UAV Dataset

The Montana VNIR dataset consists of a drone-acquired hyperspectral mosaic with 273 spectral bands and spatial dimensions of 10706 × 8360 pixels. Before applying any band-selection algorithms or training classification models, we conducted a comprehensive exploratory data analysis (EDA) to (i) construct an accurate pixel-wise ground truth (GT), (ii) assess class distribution and label sparsity, and (iii) evaluate spectral-band radiometric quality and class-wise separability.

2.4.1. Ground-Truth Construction and Alignment Verification

Field plot measurements, provided as longitude–latitude points with species labels, were converted into a pixel-wise raster mask using a custom reprojection and rasterization pipeline. The GPS coordinates (EPSG:4326) were transformed into the coordinate reference system of the VNIR cube using rasterio.warp.transform, followed by conversion to pixel indices via the GeoTIFF affine transform. Each species was assigned a unique integer class ID (with 0 reserved for background), and points within the VNIR footprint were written into a 2-D mask saved as montana_gt.npy.

To validate GT–image alignment, we generated a pseudo-RGB composite generated using three high-SNR VNIR bands (R=150, G=100, B=50) and overlaid with the GT to verify geometric alignment. This procedure revealed and resolved an initial misalignment between labels and imagery. After reconstructing the GT using accurate reprojection, the overlay confirmed correct spatial correspondence between labeled pixels and canopy/burned regions. As shown in Figure 2, the corrected GT aligns well with canopy structures visible in the hyperspectral mosaic. This verification step is essential, as even small registration errors in UAS hyperspectral imagery can significantly degrade classification accuracy.

2.4.2. Class Distribution and Label Sparsity

The final GT map contains five forest species or condition classes: Dead Biomass, Dead Ponderosa Pine, Douglas Fir, Ponderosa Pine, and Western Larch. The labels are extremely sparse relative to the full image extent, and the class distribution is highly imbalanced. Specifically, we observe 287 pixels for Dead Biomass, 212 for Ponderosa Pine, 183 for Douglas Fir, 49 for Western Larch, and only 30 for Dead Ponderosa Pine. Such an imbalance necessitates the use of class-balanced sampling, focal-loss formulations, or class weighting during model training to avoid bias toward majority classes.

2.4.3. Band-wise Radiometric Quality Analysis

Each of the 273 bands was evaluated using descriptive and noise-related statistics computed across the full spatial domain, including the mean, standard deviation, minimum, maximum, dynamic range, coefficient of variation ($\sigma /\mu$), zero-fraction (proportion of pixels equal to zero), and a simple SNR proxy ($\mu /(\sigma +ϵ)$).

The zero-fraction curve in Figure 3 reveals substantial sensor noise at the spectral extremes. Bands 1–40 show very high zero-fraction values ($0.60$–$0.68$), indicating weak short-wavelength sensor responsivity. A similar degradation is observed at the high-wavelength end (bands 240–273), where the detector suffers from end-of-range rolloff.

The standard deviation curve in Figure 4 complements this observation. Bands 40–130 exhibit stable but moderate variance, whereas bands 130–240 show the strongest spatial variability ($\approx 0.014$–$0.016$), representing the most informative spectral region for vegetation and burn discrimination. Beyond band 240, variance drops sharply again due to sensor noise.

Based on these EDA diagnostics, bands 1–40 and 240–273 were removed prior to downstream modeling, leaving 211 high-quality bands for analysis.

2.4.4. Radiometric Normalization and Clean Cube Creation

For each retained band, we applied robust percentile clipping (1st–99th percentile) followed by min–max normalization to the range $[0,1]$. This process reduces the influence of shadows and extreme outliers, while mitigating cross-track illumination variation. The resulting cleaned VNIR cube serves as the standardized input for subsequent band-selection and classification experiments.

2.4.5. Class-wise Spectral Signature Analysis

Using the cleaned cube and GT mask, we computed the mean spectral signature for each class by averaging spectral reflectance values across all labeled pixels:

$${\overline{s}}_c\left(b\right)={\displaystyle \frac{1}{N_c}}\sum _{\left\{\right(x,y):GT(x,y)=c\}}{I}_b(x,y)$$

(3)

The resulting spectra showed that the three live conifer species (Douglas Fir, Ponderosa Pine, Western Larch) share closely aligned spectral profiles, with only subtle differences. Dead biomass exhibits lower reflectance and altered mid-VNIR slopes, but separability is distributed across many bands. Unlike benchmark datasets (e.g., Indian Pines, Pavia), spectral contrast in this dataset is broad rather than concentrated, indicating that larger Top-K values (e.g., 60–100 bands) are necessary for effective band selection.

2.4.6. Montana EDA Insights

The Montana EDA demonstrates that: (i) accurate GT construction and alignment are essential for reliable supervision, (ii) extreme class imbalance requires balanced sampling strategies, (iii) only the mid-VNIR range (bands 40–240) provides high-quality spectral information, (iv) short- and long-wavelength bands should be discarded due to sensor noise, and (v) the subtle, distributed spectral differences between forest classes motivate the use of larger band subsets and patch-based deep models.

These findings guided all subsequent preprocessing, band-selection, and classification experiments.

2.5. Feature Selection (Band Selection)

To mitigate high dimensionality and spectral redundancy in hyperspectral images (HSIs), four band-selection techniques were employed: Principal Component Analysis (PCA), Spatial-Spectral Edge Preservation (SSEP), Spectral-Redundancy Penalized Attention Ranking (SRPA), and Deep Reinforcement Learning (DRL). Each method processes a pre-cleaned HSI data cube (height × width × bands, stored as .npy) and the corresponding ground truth (GT) map (2D integer labels, .npy) after noise removal and normalization. These methods include unsupervised (PCA) and label-guided or supervised (SSEP, SRPA, DRL) band selection strategies, producing ranked band indices, scores, and visualizations for subsequent classification tasks. For the Montana UAS dataset, the EDA-cleaned cube additionally includes robust percentile clipping prior to normalization, as described in Section 2.4.

2.6. Principal Component Analysis (PCA)-Based Band Selection

PCA is a classical unsupervised technique for variance-based band prioritization. We fit PCA on spectra from labeled pixels (GT > 0) to focus variance estimation on class-relevant regions. Full PCA is fitted (n_components = B, where B is the number of bands), yielding explained variance ratios (EVR) and component loadings. The score for each band b is computed as

$$\sum _{i=1}^M{\mathrm{EVR}}_i\times {\left({\mathrm{loading}}_{i,b}\right)}^2$$

(4)

We set $M=30$ as a stable trade-off that captures dominant variance while avoiding noise-dominated components; sensitivity to M was empirically minor in pilot runs [36].

2.7. Spatial-Spectral Edge Preservation (SSEP)

SSEP is a label-guided band ranking method that leverages spatial class boundaries derived from ground-truth annotations [37]. A reference binary edge map is derived from the GT using the Sobel gradient operator on label transitions. For each band, light Gaussian smoothing ($\sigma =1.0$) was applied to reduce noise, followed by Sobel edge computation, thresholding at the 95th percentile, and binarization. Alignment with the GT edge map is quantified using the Dice coefficient:

$$\mathrm{Dice}={\displaystyle \frac{2\times |E_{\mathrm{band}}\cap E_{\mathrm{GT}}|}{|E_{\mathrm{band}}|+|E_{\mathrm{GT}}|}}$$

(5)

Bands are ranked in descending order of these scores, favoring those that maintain sharp spatial-spectral edges [37]. The implementation of SSEP [32] in this study follows the methodology described in Algorithm 1.

Algorithm 1: SSEP Band Selection

Require: HSI cube C ($H\times W\times B$), GT map G ($H\times W$) Ensure: Ranked band indices order (descending scores) 1: edge_GT← binarize(Sobel(float(G))) ▹ Binary edges where gradient $>0$ 2: scores← zeros(B) 3: for$b=0$to$B-1$do 4:     band_img ←$C[:,:,b]$ 5:     smoothed ← Gaussian_filter(band_img, $\sigma =1.0$) 6:     edge_band ← Sobel(smoothed) 7:     thresh ← percentile(edge_band, 95) 8:     edge_bin ← (edge_band> thresh) 9:     scores[b]← Dice(edge_bin, edge_GT) 10: end for 11: order← argsort(scores)[::-1] 12: returnorder, scores

2.8. Spectral-Redundancy Penalized Attention Ranking (SRPA)

SRPA employs a supervised attention mechanism to balance band informativeness and diversity. Small $5\times 5$ patches are extracted around labeled pixels (up to 4000 samples per dataset; when fewer labeled samples were available, as in the Montana UAS dataset, all labeled patches were used) and split into train/validation sets (80/20, stratified). A lightweight 3D CNN (Conv3D layers, max-pooling, global average pooling, and Squeeze-and-Excitation block) is trained for two epochs to obtain stable attention trends while keeping band-ranking stage computationally lightweight. After training, mean band-wise attention weights are inferred from the SE block. Redundancy is computed from the correlation matrix of flattened subsampled patches:

$$\mathrm{red}\left(b\right)={\displaystyle \frac{1}{B-1}}\sum _{j\ne b}\left|\mathrm{corr}(b,j)\right|$$

(6)

This term penalizes bands that are highly correlated with the remainder of the spectrum. Final scores are attention $-\lambda \times$ redundancy ($\lambda =0.3$). Bands are ranked in descending order, promoting discriminative yet non-redundant selections [38]. The complete SRPA procedure, including patch extraction, network training, attention aggregation, and redundancy penalization, is summarized in Algorithm 2.

Algorithm 2: SRPA Band Selection

Require: HSI cube C ($H\times W\times B$), GT map G, patch_size=5, max_patches=4000, $\lambda =0.3$ Ensure:Ranked band indices order (descending scores) 1: Extract patches P ($N\times 5\times 5\times B$) and labels y from labeled positions ($N\le$ max_patches) 2: Split P, y→ train/validation (80/20, stratified) 3: Initialize 3D CNN with SE attention block 4: Train model on train patches (2 epochs, cross-entropy loss, Adam) 5: attn← mean(SE outputs over validation) ▹ shape (B,) 6: $X\leftarrow$ flatten(subsample(P)) ▹ ($N_{\mathrm{flat}}\times B$) 7: corr← corrcoef(X) 8: redundancy← (sum($\left|\mathtt{corr}\right|$, axis=1) − 1) / ($B-1$) 9: scores←attn$-\lambda \times$redundancy 10: order← argsort(scores)[::-1] 11: returnorder, scores, attn, redundancy

2.9. Deep Reinforcement Learning (DRL)–Based Band Selection

To enable adaptive and dataset-agnostic band selection, we formulate spectral band selection as a sequential decision-making problem and solve it using deep reinforcement learning. The objective is to incrementally construct compact band subsets that maximize downstream classification performance while discouraging redundant or unnecessarily large selections. The problem is modeled as a Markov Decision Process (MDP) and optimized using a Deep Q-Network (DQN), allowing the agent to explicitly reason about long-term selection quality rather than relying on greedy or one-shot ranking strategies.

2.9.1. Markov Decision Process Formulation

Band selection is formulated as a finite-horizon MDP defined by the tuple $(\mathcal{S},\mathcal{A},P,R,\gamma )$, where $\mathcal{S}$ denotes the state space, $\mathcal{A}$ the action space, P the transition dynamics, R the reward function, and $\gamma \in (0,1]$ the discount factor.

State.

At decision step t, the state $s_t\in \mathcal{S}$ encodes the current band-selection context and is defined as:

$$s_t=\left[{\mathbf{m}}_t,\left|S_t\right|/B,{\mathbf{r}}_t\right],$$

(7)

where $S_t$ is the set of selected bands up to step t, B denotes the number of spectral bands after EDA-based cleaning (dataset-specific; e.g., $B=211$ for the Montana VNIR dataset), ${\mathbf{m}}_t\in {\{0,1\}}^B$ is a binary mask indicating selected bands, and ${\mathbf{r}}_t\in {\mathbb{R}}^B$ represents per-band redundancy statistics computed from the EDA-cleaned hyperspectral cube. This state representation captures selection history, relative subset size, and spectral redundancy.

Action.

The action space consists of selecting one unselected spectral band:

$$a_t\in \{b\mid b\notin S_t\}.$$

(8)

Actions corresponding to previously selected bands are masked to prevent reselection.

Transition.

Transitions are deterministic. Executing action $a_t$ updates the selected set as $S_{t+1}=S_t\cup \left\{a_t\right\}$, yielding the next state $s_{t+1}$. An episode terminates when a predefined band budget K is reached.

Reward.

The reward function balances downstream classification performance and subset compactness. After selecting K bands, a lightweight classifier is trained using only the selected subset, and the terminal reward is defined as:

$$R=\mathrm{Acc}\left(S_K\right)-\lambda {\displaystyle \frac{|S_K|}{B}},$$

(9)

where $\mathrm{Acc}\left(S_K\right)$ denotes validation accuracy obtained using the selected bands and $\lambda >0$ controls subset-size regularization. Intermediate rewards are set to zero, encouraging the agent to optimize long-term selection quality.

2.9.2. Deep Q-Network Architecture

The action-value function $Q(s,a)$ is approximated using a fully connected neural network that maps the state representation to Q-values over candidate actions. The network consists of multiple dense layers with ReLU activations. A target network and an experience replay buffer are employed to stabilize training and mitigate overestimation bias.

2.9.3. Training Procedure

At each episode, the agent sequentially selects spectral bands until the budget K is reached. The terminal reward is computed using validation accuracy, and transitions $(s_t,a_t,R,s_{t+1})$ are stored in the replay buffer. The DQN is trained by minimizing the temporal-difference loss:

$$\mathcal{L}=\mathbb{E}\left[{\left(R+\gamma \underset{a^{\prime }}{max}{Q}_{\mathrm{target}}(s_{t+1},a^{\prime })-Q(s_t,a_t)\right)}^2\right].$$

(10)

An $ϵ$-greedy strategy is used to balance exploration and exploitation during training.

2.9.4. Dataset-Specific Handling

The DRL formulation is dataset-agnostic. The value of B and redundancy statistics are derived from the EDA-cleaned hyperspectral cube for each dataset. For the Montana UAV dataset, EDA cleaning includes percentile clipping and removal of low-quality VNIR bands, as described in Section 2.4. No additional dataset-specific tuning is applied.

2.9.5. Evaluation and Stability

Due to the stochastic nature of reinforcement learning, each DRL experiment is repeated multiple times with different random seeds. Band rankings are obtained by averaging selection frequencies across runs, yielding stable and reproducible band importance estimates.

The DRL-based band selection method enables adaptive, sequential selection of informative spectral bands while explicitly accounting for subset compactness. By integrating EDA-derived statistics and downstream classification feedback into a unified MDP framework, the method provides a flexible alternative to classical and attention-based band-ranking strategies.

2.10. Classification Models

To evaluate the effectiveness of the selected hyperspectral bands, we employed a suite of classical machine learning classifiers and a deep learning model. These models were trained on the reduced-dimensionality data obtained from the band selection methods (PCA, SSEP, SRPA, and DRL) as well as the full-band baseline. The classifiers include Random Forest (RF), K-Nearest Neighbors (KNN), Support Vector Machines (SVM), and a 3D Convolutional Neural Network (3D-CNN). These models were chosen for their proven efficacy in hyperspectral image classification tasks, balancing computational efficiency with performance in handling high-dimensional spectral-spatial data [35,39]. All models were implemented using scikit-learn for classical methods and PyTorch for the deep model, with hyperparameters tuned based on standard practices in HSI literature [40].

2.10.1. Random Forest

RF is an ensemble learning method that constructs multiple decision trees during training and outputs the class that is the mode of the classes from individual trees. In our implementation, we used the RandomForest Classifier from scikit-learn with 200 estimators, utilizing all available CPU cores (n_jobs=-1) for parallel processing, and a fixed random state of 0 for reproducibility. RF is particularly effective for HSI classification due to its ability to handle multicollinearity in spectral bands, and its robustness to overfitting [41]. The model was trained on flattened spectral features from labeled pixels or selected bands, without incorporating spatial information explicitly. This approach aligns with prior studies on hyperspectral data where RF serves as a strong baseline for pixel-wise classification [42].

2.10.2. KNN

The KNN classifier assigns a class to a query point based on the majority vote of its k nearest neighbors in the feature space. We implemented KNN using scikit-learn’s KNeighborsClassifier with k=5, employing the default Euclidean distance metric. This value of k was selected based on empirical performance in HSI datasets, where smaller neighborhoods help capture local spectral similarities while avoiding excessive noise sensitivity [43]. KNN is computationally simple and non-parametric, making it suitable for hyperspectral data with varying class distributions, as seen in post-fire fuel classification scenarios [44]. Like RF, it operates on pixel-based spectral features, leveraging the reduced band subsets to mitigate the curse of dimensionality.

2.10.3. SVM

SVM aims to find the optimal hyperplane that separates classes in a high-dimensional space, maximizing the margin between support vectors. Our SVM implementation used scikit-learn’s SVC with a radial basis function (RBF) kernel, regularization parameter C=10, and gamma set to ’scale’ (automatically computed as 1/($n_{f}eatures$ * X.var())). This configuration is commonly used in HSI classification to handle non-linear separability in spectral data [45]. SVM’s strength lies in its effectiveness with small sample sizes and high-dimensional inputs, which is relevant for our prescribed-burn dataset where labeled pixels may be sparse [46]. The model was trained on the same pixel-wise features as RF and KNN, benefiting from band selection to reduce kernel computation overhead.

2.10.4. 3D-CNN

For all datasets, spatial–spectral patches of fixed size were extracted from the EDA-cleaned hyperspectral cube prior to 3D-CNN training, ensuring consistent spatial context across band-selection methods. The 3D Convolutional Neural Network (3D-CNN) extends traditional CNNs by incorporating spectral depth as an additional dimension, enabling joint spatial-spectral feature extraction. Our model, implemented in PyTorch, consists of two 3D convolutional layers followed by batch normalization, ReLU activation, max pooling, adaptive average pooling, and a fully connected output layer. The input shape is ($batc{h}_{s}ize$, 1, bands, $patc{h}_{s}ize$, $patc{h}_{s}ize$), where patches are centered on labeled pixels with a $patc{h}_{s}ize$ of 5 (as defined in the dataset configurations). The first convolutional layer uses 16 filters with a kernel size of (3,3,7) and padding to preserve dimensions, while the second uses 32 filters with (3,3,5). Training was performed over 20 epochs with Adam optimization (learning rate 1e-3), cross-entropy loss, and a batch size of 32, monitoring validation accuracy to select the best model [47]. This architecture is tailored for HSI, capturing volumetric patterns in post-fire environments, and has shown superior performance over 2D CNNs in spectral-spatial tasks [48].

2.11. Evaluation Metrics

Model performance was assessed using a comprehensive set of metrics to ensure balanced evaluation across imbalanced classes typical in HSI datasets. Overall Accuracy (OA) measures the percentage of correctly classified pixels. Cohen’s Kappa coefficient accounts for chance agreement, providing a more robust measure of inter-class reliability [49]. Macro-averaged Precision, Recall, and F1-score were computed to evaluate performance across all classes equally, handling class imbalances in prescribed-burn scenarios where certain fuel types (e.g., char or recovering vegetation) may be underrepresented [50]. These metrics were calculated in percentage form for consistency with OA, using scikit-learn’s $precisio{n}_{s}core$, $recal{l}_{s}core$, and $f{1}_{s}core$ functions with ’macro’ averaging and $zer{o}_{d}ivision$=0 to manage undefined cases. Additionally, confusion matrices were generated to visualize misclassifications. For deep models, the best validation accuracy during training was also logged as an indicator of convergence [51]. All metrics were derived from a 70/30 train-validation split, stratified by class to maintain distribution. All classifier hyperparameters were fixed to commonly used values from prior hyperspectral classification literature and were kept identical across all band-selection methods to ensure a fair comparative evaluation. All reported metrics were computed exclusively on held-out validation or test data and were not used during training or band-selection stages.

2.12. Experimental Setup

Experiments were conducted on a unified pipeline implemented in Python 3.12, utilizing NumPy for data handling, scikit-learn for classical models, and PyTorch for deep architectures. Datasets included benchmark HSI scenes (Indian Pines, Pavia University, Salinas, Botswana, KSC) and a custom VNIR hyperspectral dataset from prescribed burns at Lubrecht Experimental Forest, Montana, USA. Data preprocessing involved loading hyperspectral cubes and ground truth maps as NumPy arrays, extracting labeled pixels or center patches ($patc{h}_{s}ize$=5 for patch-based models), and applying a 70/30 train-validation split with stratification ($rando{m}_{s}tate$=0). For band selection, methods (PCA, SSEP, SRPA, DRL) were applied to generate top-k band subsets (k in [5,10,15,20,25,30,35,40,45,50]), with results compared against a full-band baseline. Training was performed on an NVIDIA GPU-enabled system where available, with CPU fallback for classical models. Results were logged to a CSV file tracking dataset, method, top-k, model, and metrics. Cross-dataset evaluation ensured generalizability, particularly testing Montana data for real-world fire monitoring applicability [52]. The setup emphasizes reproducibility, with all paths and hyperparameters defined in a central configuration file. All experiments were repeated over multiple random seeds, and results were averaged to reduce variance.

In addition to dataset-specific evaluations, we report a cross-dataset summary in which, for each band selection technique, the best-performing configuration is selected based on aggregated performance across all datasets. This analysis is intended to highlight general performance trends and preprocessing sensitivity rather than dataset-specific accuracy.

3. Results

3.1. Result Aggregation and Best-Configuration Identification

For each dataset, all evaluated configurations—defined by the combination of band selection technique (No Band Selection, PCA, SSEP, SRPA, DRL), number of selected bands (Top-K), and classifier (RF, SVM, KNN, 3DCNN)—were executed using five independent runs with different random seeds. Performance is reported as mean ± standard deviation for Overall Accuracy (OA), Cohen’s Kappa, and macro-averaged F1 score.

This exhaustive evaluation produced approximately 400–415 configuration-level results per dataset, which are provided as Supplementary Material. To present a concise yet representative summary in the main manuscript, we selected one best-performing configuration per band selection technique under EDA and No-EDA conditions, respectively.

Best configurations were identified using a hierarchical selection criterion that prioritizes: (i) highest mean OA, (ii) highest mean Kappa, (iii) highest mean macro-F1, and (iv) lower OA standard deviation. The resulting best-per-technique summaries are reported in Tables 1–6, with Table 7 providing a consolidated cross-dataset comparison.

3.2. Dataset-Specific Results

3.2.1. Indian Pines (Table 2)

Indian Pines exhibits strong sensitivity to both band selection strategy and preprocessing. As shown in Table 2, DRL achieves the highest OA, Kappa, and macro-F1 under both EDA and No-EDA conditions, typically selecting compact Top-K band subsets combined with 3D-CNN classifiers. SRPA consistently ranks second, also favoring reduced band subsets, while PCA yields the lowest performance across all metrics.

The $\Delta$ columns in Table 2 show positive $\Delta$OA and $\Delta$F1 values for DRL, SRPA, and SSEP, indicating that EDA improves the best-achievable performance of learning-based techniques on this spectrally complex and class-imbalanced dataset.

3.2.2. Pavia University (Table 3)

Pavia University shows a clear dominance of DRL, which achieves the highest OA, Kappa, and macro-F1 under both EDA and No-EDA conditions (Table 3). Best-performing configurations consistently involve small Top-K subsets and CNN-based classifiers, highlighting the importance of spatial–spectral modeling in urban environments.

EDA provides consistent but moderate gains, as reflected by positive $\Delta$OA and $\Delta$F1 values for DRL and SRPA. PCA again exhibits substantially lower performance compared to learning-based methods.

3.2.3. Salinas (Table 4)

Salinas achieves high OA across all techniques (Table 4). DRL produces the best performance, followed closely by SRPA, with best configurations typically involving reduced band subsets and RF or CNN-based classifiers.

The $\Delta$ values in Table 4 are small but predominantly positive, suggesting that EDA offers incremental improvements in best-achievable performance without altering the overall ranking of techniques.

3.2.4. Botswana (Table 5)

Botswana demonstrates uniformly high classification performance across all techniques (Table 5). DRL marginally outperforms other methods, with best configurations typically involving moderate Top-K values and classical classifiers such as SVM or RF. Performance differences between techniques are small.

The $\Delta$ values in Table 5 are close to zero for all techniques, indicating that EDA has minimal impact on the best-achievable performance. This behavior is consistent with the dataset’s relatively clean spectral signatures and strong class separability.

3.2.5. Kennedy Space Center (KSC) (Table 6)

For the KSC dataset, DRL achieves the highest OA and macro-F1, followed by SRPA, while PCA again performs worst (Table 6). Best configurations frequently involve intermediate Top-K values combined with RF or SVM classifiers.

The $\Delta$ values in Table 6 are mixed across techniques, with some methods showing small positive differences and others exhibiting near-zero or negative values. This indicates that the impact of EDA on KSC is technique-dependent rather than uniformly beneficial.

3.2.6. Montana (UAV VNIR) (Table 7)

Montana exhibits fundamentally different behavior compared to benchmark datasets. As shown in Table 7, overall OA and macro-F1 are substantially lower, and standard deviations are higher, indicating reduced stability across runs.

In contrast to other datasets, baseline or classical methods remain competitive, while DRL and SRPA do not consistently dominate. Best configurations often retain larger Top-K values and favor classical classifiers. The $\Delta$ values are inconsistent in both sign and magnitude, indicating that EDA does not reliably improve best-achievable performance for this real-world UAV VNIR dataset.

3.3. Cross-Dataset Trends (Table 8)

While dataset-specific results provide insight into method performance under individual data characteristics, it is also informative to examine trends that hold consistently across datasets. To this end, we analyze cross-dataset behavior by aggregating results across all benchmark datasets and selecting, for each band selection technique, the globally best-performing configuration based on overall classification accuracy.

Table 8 summarizes the impact of Exploratory Data Analysis (EDA) under these globally selected configurations. Although absolute performance varies across datasets, the aggregated results reveal consistent method-level trends. In particular, DRL-based band selection benefits most from EDA, exhibiting substantial improvements in overall accuracy, Kappa coefficient, and macro-F1 score. PCA-based selection shows moderate but consistent gains, indicating sensitivity to preprocessing and normalization. In contrast, edge-based methods such as SSEP and SRPA demonstrate mixed behavior, with marginal or negative performance changes under certain configurations.

These observations suggest that EDA plays a more critical role for learning-based and sequential band selection methods, which are more sensitive to spectral noise, scale variation, and redundancy. The cross-dataset analysis complements the dataset-specific evaluations and supports the robustness of the observed trends across heterogeneous hyperspectral scenes.

3.4. Summary of Results

Overall, the results in Tables 2–7 demonstrate that learning-based band selection methods, particularly DRL, consistently achieve the strongest performance on benchmark hyperspectral datasets, while EDA provides dataset-dependent improvements. The Montana dataset further illustrates that these gains do not universally transfer to real-world conditions, underscoring the importance of dataset characteristics when selecting band- selection and preprocessing strategies.

4. Discussion

The results presented in this study highlight several important insights regarding hyperspectral band-selection, preprocessing strategies, and model behavior across both benchmark and real-world datasets. Rather than focusing solely on absolute performance, this discussion interprets the observed trends in the context of dataset characteristics, model complexity, and the role of exploratory data analysis (EDA).

4.1. Effectiveness of Learning-Based Band Selection

Across five of the six evaluated datasets, learning-based band selection methods—particularly deep reinforcement learning (DRL)—achieved the highest classification performance. This consistent dominance suggests that DRL is effective at jointly optimizing spectral relevance and downstream classification performance, especially when sufficient labeled data and relatively stable spectral signatures are available. Unlike classical techniques such as PCA, which prioritize variance preservation, DRL explicitly optimizes task-specific objectives, leading to improved class discrimination as reflected by higher macro-F1 and Kappa scores.

The attention-based SRPA method consistently ranked second across all datasets, demonstrating a strong balance between performance and stability. Compared to DRL, SRPA appears less sensitive to dataset size and noise, which may explain its competitive behavior even in more challenging settings. These observations support the growing body of literature advocating task-aware and learning-based band-selection methods over purely statistical approaches.

4.2. Limitations of Variance-Based Selection

PCA yielded the lowest performance across all datasets, confirming that variance maximization alone is insufficient for hyperspectral classification. High-variance bands do not necessarily correspond to class-discriminative features, particularly in datasets with class imbalance or overlapping spectral signatures. The consistently poor performance of PCA across both EDA and No-EDA settings reinforces the need for supervised or task-driven band selection strategies in practical hyperspectral applications.

4.3. Role of Exploratory Data Analysis

The impact of EDA was found to be strongly dataset-dependent. On spectrally complex benchmark datasets such as Indian Pines and Pavia University, EDA consistently improved best-achievable performance for learning-based methods. These improvements likely arise from the removal of noisy samples, correction of spectral outliers, and improved feature normalization, which collectively enhance model generalization and stability across runs.

In contrast, datasets with cleaner spectral characteristics, such as Botswana and Salinas, showed minimal performance differences between EDA and No-EDA. Importantly, EDA did not degrade performance on any dataset, indicating that it is generally a safe preprocessing step, though its benefits may be marginal when data quality is already high.

4.4. Challenges of Real-World UAV Hyperspectral Data

The Montana UAV VNIR dataset exhibited fundamentally different behavior compared to benchmark datasets. Overall classification performance was substantially lower, and learning-based methods did not consistently outperform simpler baseline or classical models. This divergence can be attributed to several real-world challenges, including limited labeled samples, increased spectral noise, class imbalance, and domain-specific variability introduced by UAV acquisition conditions.

In this setting, EDA did not provide consistent performance gains, suggesting that preprocessing alone cannot compensate for insufficient supervision or highly variable data distributions. These findings underscore an important limitation of deep and reinforcement-based band selection methods: while powerful, they remain sensitive to data quality and training signal strength.

4.5. Implications for Practical Deployment

Taken together, the results suggest that the choice of band selection strategy should be guided by dataset characteristics rather than assumed universally optimal. Learning-based methods such as DRL and SRPA are highly effective on well-annotated benchmark datasets and structured scenes, but simpler methods may remain competitive in real-world scenarios with limited labels. The dataset-dependent impact of EDA further emphasizes the need for adaptive preprocessing pipelines rather than fixed workflows.

4.6. Future Research Directions

Future work will focus on improving the robustness of learning-based band selection methods for real-world UAV hyperspectral data. Promising directions include incorporating label-efficient or self-supervised learning strategies, leveraging spatial context more effectively, and designing hybrid approaches that combine the stability of classical methods with the adaptability of learning-based models. Extending the evaluation framework to additional sensing modalities and ecological monitoring tasks will further support the generalization of these findings.

5. Conclusions

This work presented a comprehensive and reproducible evaluation of hyperspectral band selection techniques across six datasets, including widely used benchmark scenes and a real-world UAV-based VNIR dataset. Consistent with the objectives stated in the Abstract, we systematically compared classical, attention-based, and reinforcement learning–based band selection methods under both EDA and No-EDA settings using a unified experimental protocol.

By aggregating results over five independent runs and selecting representative configurations based on mean performance and stability, we ensured a fair and transparent comparison across techniques. The results demonstrate that learning-based band selection methods consistently outperform classical approaches on benchmark datasets. In particular, deep reinforcement learning (DRL) achieved the highest overall accuracy on five out of six datasets, confirming its robustness across diverse hyperspectral scenes. The attention-based SRPA method also showed strong and stable performance, frequently ranking second and consistently outperforming PCA-based selection.

The effect of exploratory data analysis (EDA) was found to be strongly dataset-dependent. On spectrally complex benchmarks such as Indian Pines and Pavia University, EDA improved the best-achievable performance for learning-based methods, as reflected by positive $\Delta$ values in overall accuracy and macro-F1. On cleaner datasets such as Botswana and Salinas, the impact of EDA was marginal but non-degrading. These findings support the Abstract’s claim that preprocessing can enhance performance, but should not be assumed to be universally beneficial.

In contrast to benchmark datasets, the Montana UAV VNIR dataset exhibited substantially lower performance across all methods, with higher variability across runs. In this real-world setting, learning-based band selection methods did not consistently outperform simpler baseline or classical approaches, and EDA did not provide reliable improvements. This observation highlights the limitations of deep and reinforcement-based methods when applied to noisy data with limited labeled samples, and directly supports the Abstract’s emphasis on the challenges of real-world hyperspectral analysis.

Overall, this study confirms that learning-based band selection—particularly DRL—offers significant advantages for hyperspectral image classification when data quality and supervision are sufficient. At the same time, the observed variability across datasets underscores the importance of aligning band selection strategies, preprocessing choices, and model complexity with dataset characteristics. Future work will focus on improving robustness for real-world UAV hyperspectral data, incorporating label-efficient and self-supervised learning strategies, and extending the proposed evaluation framework to additional sensing modalities and ecological applications.