# Augmenting Heritage Ocean-Color Aerosol Models for Enhanced Remote Sensing of Inland and Nearshore Coastal Waters

Martin Montes, Nima Pahlevan, David Giles, Jean-Claude Roger, Peng-Wang Zhai, Brandon Smith, Robert Levy, P. Jeremy Werdell, Alexander Smirnov

► To cite this version:

Martin Montes, Nima Pahlevan, David Giles, Jean-Claude Roger, Peng-Wang Zhai, et al. Augmenting Heritage Ocean-Color Aerosol Models for Enhanced Remote Sensing of Inland and Nearshore Coastal Waters. Frontiers in Remote Sensing, 2022, 3, 10.3389/frsen.2022.860816. hal-04515795

HAL Id: hal-04515795

https://uca.hal.science/hal-04515795

Submitted on 16 May 2024

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>
<figure><img src="image_1.png" /></figure>

Check for
updates

# Augmenting Heritage Ocean-Color Aerosol Models for Enhanced Remote Sensing of Inland and Nearshore Coastal Waters

Martin Montes<fnref n="1" />,², Nima Pahlevan<fnref n="1" />,²*, David M. Giles<fnref n="1" />,², Jean-Claude Roger²,³, Peng-wang Zhai<fnref n="4" />, Brandon Smith<fnref n="1" />,², Robert Levy², P. Jeremy Werdell² and Alexander Smirnov<fnref n="1" />,²

Satellite remote sensing of near-surface water composition in terrestrial and coastal regions is challenging largely due to uncertainties linked to a lack of representative continental aerosols in the atmospheric correction (AC) framework. A comprehensive family of absorbing aerosols is proposed by analyzing global AERONET measurements using the Partition Around Medoids (PAM) classifier. The input to the classifier is composed of Version 3, Level 2.0 daily average aerosol properties [i.e., single scattering albedo at $λ = 0.44$ µm, (SSA(0.44)) and the Angstrom exponents for extinction and absorption $AE_e(0.44–0.87)$ and $AE_a(0.44–0.87)$, respectively from observations from June 1993 to September 2019. The PAM classification based on low daily aerosol optical depth ($AOD(0.44) ≤ 0.4$) suggested 27 distinct aerosol clusters encompassing five major absorbing aerosol types (Dust (DU), Marine (MAR), Mixed (MIX), Urban/Industrial (U/I), and Biomass Burning (BB)). Seasonal patterns of dominant PAM-derived clusters at three AERONET sites (GSFC, Kanpur, and Banizoumbou) strongly influenced by U/I, DU, and BB types, respectively, showed a satisfactory agreement with variations of aerosol mixtures reported in the literature. These PAM-derived models augment the National Aeronautics and Space Administration's (NASA's) aerosol models (A2010) applied in its operational AC. To demonstrate the validity and complementary nature of our models, a coupled ocean-atmosphere radiative transfer code is employed to create a simulated dataset for developing two experimental machine-learning AC processors. These two processors differ only in their aerosol models used in training: 1) a processor trained with the A2010 aerosol models (ACI) and 2) a processor trained with both PAM and A2010 aerosol models (ACII). These processors are applied to Landsat-8 Operational Land Imager (OLI) matchups ($N = 173$) from selected AERONET sites equipped with ocean color radiometers (AERONET-OC). Our assessments showed improvements of up to 30% in retrieving remote sensing reflectance ($R_{rs}$) in the blue bands. In general, our empirically derived PAM aerosol models complement A2010 models (designed for regions strongly influenced by marine conditions) over continental and coastal waters where absorbing aerosols are present (e.g., urban environments, areas impacted by dust, or wildfire events).

**OPEN ACCESS**

**Edited by:**
Jianwei Wei,
Center for Satellite Applications and
Research (STAR), United States

**Reviewed by:**
Lei Bi,
Zhejiang University, China
Jae-Hyun Ahn,
Korea Institute of Ocean Science and
Technology (KIOST), South Korea

*Correspondence:
Nima Pahlevan
nima.pahlevan@nasa.gov

**Specialty section:**
This article was submitted to
Multi- and Hyper-Spectral Imaging,
a section of the journal
Frontiers in Remote Sensing

Received: 23 January 2022

Accepted: 16 March 2022

Published: 02 May 2022

**Citation:**
Montes M, Pahlevan N, Giles DM, Roger J-C, Zhai P-w, Smith B, Levy R, Werdell PJ and Smirnov A (2022) Augmenting Heritage Ocean-Color Aerosol Models for Enhanced Remote Sensing of Inland and Nearshore Coastal Waters. *Front. Remote Sens.* 3:860816
doi: 10.3389/frsen.2022.860816

<footnote n="1">¹Science Systems and Applications, Inc., (SSAI), Lanham, MD, United States, ²NASA Goddard Space Flight Center, Greenbelt, MD, United States, ³Department of Geographical Sciences, University of Maryland, College Park, MD, United States,</footnote>

<footnote n="4">⁴Department of Physics, University of Maryland Baltimore County, Baltimore, MD, United States</footnote>
<figure><img src="image_2.png" /></figure>

Check for
updates

# Augmenting Heritage Ocean-Color Aerosol Models for Enhanced Remote Sensing of Inland and Nearshore Coastal Waters

Martin Montes<fnref n="1" />,², Nima Pahlevan<fnref n="1" />,²*, David M. Giles<fnref n="1" />,², Jean-Claude Roger²,³, Peng-wang Zhai<fnref n="4" />, Brandon Smith<fnref n="1" />,², Robert Levy², P. Jeremy Werdell² and Alexander Smirnov<fnref n="1" />,²

Satellite remote sensing of near-surface water composition in terrestrial and coastal regions is challenging largely due to uncertainties linked to a lack of representative continental aerosols in the atmospheric correction (AC) framework. A comprehensive family of absorbing aerosols is proposed by analyzing global AERONET measurements using the Partition Around Medoids (PAM) classifier. The input to the classifier is composed of Version 3, Level 2.0 daily average aerosol properties [i.e., single scattering albedo at $λ = 0.44$ µm, (SSA(0.44)) and the Angstrom exponents for extinction and absorption $AE_e(0.44–0.87)$ and $AE_a(0.44–0.87)$, respectively from observations from June 1993 to September 2019. The PAM classification based on low daily aerosol optical depth ($AOD(0.44) ≤ 0.4$) suggested 27 distinct aerosol clusters encompassing five major absorbing aerosol types (Dust (DU), Marine (MAR), Mixed (MIX), Urban/Industrial (U/I), and Biomass Burning (BB)). Seasonal patterns of dominant PAM-derived clusters at three AERONET sites (GSFC, Kanpur, and Banizoumbou) strongly influenced by U/I, DU, and BB types, respectively, showed a satisfactory agreement with variations of aerosol mixtures reported in the literature. These PAM-derived models augment the National Aeronautics and Space Administration’s (NASA’s) aerosol models (A2010) applied in its operational AC. To demonstrate the validity and complementary nature of our models, a coupled ocean-atmosphere radiative transfer code is employed to create a simulated dataset for developing two experimental machine-learning AC processors. These two processors differ only in their aerosol models used in training: 1) a processor trained with the A2010 aerosol models (ACI) and 2) a processor trained with both PAM and A2010 aerosol models (ACII). These processors are applied to Landsat-8 Operational Land Imager (OLI) matchups ($N = 173$) from selected AERONET sites equipped with ocean color radiometers (AERONET-OC). Our assessments showed improvements of up to 30% in retrieving remote sensing reflectance ($R_{rs}$) in the blue bands. In general, our empirically derived PAM aerosol models complement A2010 models (designed for regions strongly influenced by marine conditions) over continental and coastal waters where absorbing aerosols are present (e.g., urban environments, areas impacted by dust, or wildfire events).

**OPEN ACCESS**

**Edited by:**
Jianwei Wei,
Center for Satellite Applications and
Research (STAR), United States

**Reviewed by:**
Lei Bi,
Zhejiang University, China
Jae-Hyun Ahn,
Korea Institute of Ocean Science and
Technology (KIOST), South Korea

*Correspondence:
Nima Pahlevan
nima.pahlevan@nasa.gov

**Specialty section:**
This article was submitted to
Multi- and Hyper-Spectral Imaging,
a section of the journal
Frontiers in Remote Sensing

Received: 23 January 2022

Accepted: 16 March 2022

Published: 02 May 2022

**Citation:**
Montes M, Pahlevan N, Giles DM, Roger J-C, Zhai P-w, Smith B, Levy R, Werdell PJ and Smirnov A (2022) Augmenting Heritage Ocean-Color Aerosol Models for Enhanced Remote Sensing of Inland and Nearshore Coastal Waters. *Front. Remote Sens.* 3:860816
doi: 10.3389/frsen.2022.860816

<footnote n="1">¹Science Systems and Applications, Inc., (SSAI), Lanham, MD, United States, ²NASA Goddard Space Flight Center, Greenbelt, MD, United States, ³Department of Geographical Sciences, University of Maryland, College Park, MD, United States,</footnote>

<footnote n="4">⁴Department of Physics, University of Maryland Baltimore County, Baltimore, MD, United States</footnote>
With the expected geographic expansion of *in situ* aquatic validation networks (e.g., AERONET-OC), the advantages of our models will be accentuated, particularly in the ultraviolet and short blue bands.

**Keywords: absorbing aerosols, aquatic remote sensing, water quality, atmospheric correction, lakes, rivers, coastal waters**

## INTRODUCTION

The application of moderate-to-high resolution imagery from satellite platforms (e.g., Landsat-8 and Sentinel-2) for studying inland water bodies (e.g., lakes, reservoirs, rivers) and nearshore coastal ecosystems has become more prevalent over the past decade (Alikas et al., 2015; Barnes et al., 2015; Gernez et al., 2017; Snyder et al., 2017; Maciel et al., 2019; Page et al., 2019; Balasubramanian et al., 2020; Hakimdavar et al., 2020; Cao and Tzortziou, 2021). This increased usage has driven the aquatic remote sensing community to formulate and implement suitable processing techniques and algorithms to generate reliable in-water products for a host of scientific and monitoring exercises (Franz et al., 2015; Watanabe et al., 2015; Saulquin et al., 2016; Slonecker et al., 2016; De Keukelaere et al., 2018; Steinmetz and Ramon, 2018; Poddar et al., 2019; Vanhellemont, 2019; Balasubramanian et al., 2020; Manuel et al., 2020; Fan et al., 2021).

One of the main sources of uncertainties in Landsat-8/Sentinel-2 data processing is the atmospheric correction (AC) (IOCCG, 2019; Pahlevan et al., 2021a) for which inaccurate representation of aerosol radiative properties is one major weakness, particularly when absorbing aerosols are present (IOCCG, 2010; Pahlevan et al., 2017a; Frouin et al., 2019). This is explained by the contribution of coupled aerosol scattering and absorption and their high spatial and temporal variability (Prospero et al., 1983). Given a particular location and time, realistic aerosol mixtures are critical for executing a robust AC on aquatic satellite scenes to retrieve high-quality spectral remote sensing reflectance ($R_{rs}$) products (Mobley, 1999), from which bio-geophysical products are derived. In particular, reducing uncertainties in remotely sensed biogeochemical variables (e.g., chlorophyll-a concentration; Chla (Pahlevan et al., 2021a)) bears significant implications on improving hydrological/hydrodynamic models (Thewes et al., 2020; Pahlevan et al., 2012), large-scale studies (Ho et al., 2019; Adrian et al., 2009), and water-related monitoring practices (IOCCG, 2018). Current estimates of uncertainties in satellite-derived $R_{rs}$ products range approximately from 15 to > 60% depending on the spectral band and AC processor (Figure 7 in Pahlevan et al. (2021a)). These uncertainties, commonly larger in the blue bands especially when assessing $R_{rs}$ over inland waters, preclude reliable retrievals of biogeochemical properties and inherent optical properties that are pivotal to quantifying desired higher-level products, such as particle composition and phytoplankton properties. Recent studies have also demonstrated the utility of the ultraviolet (UV) bands for phytoplankton characterizations and biomass estimations off the coasts of California (Kahru et al., 2021), a spectral region where effects

of absorbing aerosols are most confounding (IOCCG, 2010). The overarching purpose of this contribution is to reduce one major uncertainty component, among others [e.g., skyglint (Gilerson et al., 2018)], in the AC framework and enhance the quality of retrieved $R_{rs}$ products, particularly in the blue bands (Gordon et al., 1997; IOCCG, 2010), over inland and coastal waters.

To date, there is no real-time integration of independent aerosol datasets into AC frameworks. This limitation is due to the lack of concurrent *in situ* measurements of aerosol properties and/or poor spatial and temporal representation of atmospheric models. Indeed, the data coverage by terrestrial observation networks (e.g., the Aerosol Robotic Network; AERONET) is generally sparse despite a significant increase in the number of sampling sites worldwide during the last two decades. Likewise, simulated aerosol distributions (e.g., The Modern-Era Retrospective analysis for Research and Applications; MERRA-2) (Gelaro et al., 2017) are unable to resolve the composition of aerosol mixtures at scales comparable to those used by moderate-resolution spaceborne imagers (footprint ~ 30 m), such as the Operational Land Imager (OLI) aboard Landsat-8. Due to these caveats, NASA's approach to the characterization of aerosol properties for each satellite observation is currently carried out based on discrete aerosol categories selected dynamically (e.g., as a function of relative humidity (RH), Ahmad et al., 2010, hereafter A2010) using observed near- and shortwave-infrared reflectances (hereafter NIR and SWIR, respectively). Other relevant AC processors apply different techniques or aerosol models to account for the aerosol contribution (Pahlevan et al., 2021a), yet, we primarily focus on NASA's AC procedure adopted by the U.S. Geological Survey for Landsat-8 data processing (Franz et al., 2015; Pahlevan et al., 2017b).

The A2010 model developments are guided experimentally via observations at selected AERONET sites and are categorized according to relative humidity (RH) (Ahmad et al., 2010). The interplay between RH and aerosol optical properties has long been recognized in marine environments (Shettle and Fenn, 1979; Gerber, 1985). However, these relationships may not apply to inland and nearshore coastal waters having a weak marine influence due to the major impact of land processes (e.g., evapotranspiration (Raoufi and Beighley, 2017)) on RH. Unlike empirical models, theoretical A2010 models use a suite of AERONET measurements collected across several small open-ocean islands and three stations along the Chesapeake Bay region to constrain the aerosol properties in marine environments suited for heritage ocean-color observations (open water and coastal, respectively). Further, A2010 models have two main characteristics: 1) they have been developed by assuming two size distribution modes (i.e., fine and coarse) having an extreme range of contributions in terms of number density, and 2) they
With the expected geographic expansion of *in situ* aquatic validation networks (e.g., AERONET-OC), the advantages of our models will be accentuated, particularly in the ultraviolet and short blue bands.

**Keywords: absorbing aerosols, aquatic remote sensing, water quality, atmospheric correction, lakes, rivers, coastal waters**

## INTRODUCTION

The application of moderate-to-high resolution imagery from satellite platforms (e.g., Landsat-8 and Sentinel-2) for studying inland water bodies (e.g., lakes, reservoirs, rivers) and nearshore coastal ecosystems has become more prevalent over the past decade (Alikas et al., 2015; Barnes et al., 2015; Gernez et al., 2017; Snyder et al., 2017; Maciel et al., 2019; Page et al., 2019; Balasubramanian et al., 2020; Hakimdavar et al., 2020; Cao and Tzortziou, 2021). This increased usage has driven the aquatic remote sensing community to formulate and implement suitable processing techniques and algorithms to generate reliable in-water products for a host of scientific and monitoring exercises (Franz et al., 2015; Watanabe et al., 2015; Saulquin et al., 2016; Slonecker et al., 2016; De Keukelaere et al., 2018; Steinmetz and Ramon, 2018; Poddar et al., 2019; Vanhellemont, 2019; Balasubramanian et al., 2020; Manuel et al., 2020; Fan et al., 2021).

One of the main sources of uncertainties in Landsat-8/Sentinel-2 data processing is the atmospheric correction (AC) (IOCCG, 2019; Pahlevan et al., 2021a) for which inaccurate representation of aerosol radiative properties is one major weakness, particularly when absorbing aerosols are present (IOCCG, 2010; Pahlevan et al., 2017a; Frouin et al., 2019). This is explained by the contribution of coupled aerosol scattering and absorption and their high spatial and temporal variability (Prospero et al., 1983). Given a particular location and time, realistic aerosol mixtures are critical for executing a robust AC on aquatic satellite scenes to retrieve high-quality spectral remote sensing reflectance ($R_{rs}$) products (Mobley, 1999), from which bio-geophysical products are derived. In particular, reducing uncertainties in remotely sensed biogeochemical variables (e.g., chlorophyll-a concentration; Chla (Pahlevan et al., 2021a)) bears significant implications on improving hydrological/hydrodynamic models (Thewes et al., 2020; Pahlevan et al., 2012), large-scale studies (Ho et al., 2019; Adrian et al., 2009), and water-related monitoring practices (IOCCG, 2018). Current estimates of uncertainties in satellite-derived $R_{rs}$ products range approximately from 15 to > 60% depending on the spectral band and AC processor (Figure 7 in Pahlevan et al. (2021a)). These uncertainties, commonly larger in the blue bands especially when assessing $R_{rs}$ over inland waters, preclude reliable retrievals of biogeochemical properties and inherent optical properties that are pivotal to quantifying desired higher-level products, such as particle composition and phytoplankton properties. Recent studies have also demonstrated the utility of the ultraviolet (UV) bands for phytoplankton characterizations and biomass estimations off the coasts of California (Kahru et al., 2021), a spectral region where effects

of absorbing aerosols are most confounding (IOCCG, 2010). The overarching purpose of this contribution is to reduce one major uncertainty component, among others [e.g., skyglint (Gilerson et al., 2018)], in the AC framework and enhance the quality of retrieved $R_{rs}$ products, particularly in the blue bands (Gordon et al., 1997; IOCCG, 2010), over inland and coastal waters.

To date, there is no real-time integration of independent aerosol datasets into AC frameworks. This limitation is due to the lack of concurrent *in situ* measurements of aerosol properties and/or poor spatial and temporal representation of atmospheric models. Indeed, the data coverage by terrestrial observation networks (e.g., the Aerosol Robotic Network; AERONET) is generally sparse despite a significant increase in the number of sampling sites worldwide during the last two decades. Likewise, simulated aerosol distributions (e.g., The Modern-Era Retrospective analysis for Research and Applications; MERRA-2) (Gelaro et al., 2017) are unable to resolve the composition of aerosol mixtures at scales comparable to those used by moderate-resolution spaceborne imagers (footprint ~ 30 m), such as the Operational Land Imager (OLI) aboard Landsat-8. Due to these caveats, NASA's approach to the characterization of aerosol properties for each satellite observation is currently carried out based on discrete aerosol categories selected dynamically (e.g., as a function of relative humidity (RH), Ahmad et al., 2010, hereafter A2010) using observed near- and shortwave-infrared reflectances (hereafter NIR and SWIR, respectively). Other relevant AC processors apply different techniques or aerosol models to account for the aerosol contribution (Pahlevan et al., 2021a), yet, we primarily focus on NASA's AC procedure adopted by the U.S. Geological Survey for Landsat-8 data processing (Franz et al., 2015; Pahlevan et al., 2017b).

The A2010 model developments are guided experimentally via observations at selected AERONET sites and are categorized according to relative humidity (RH) (Ahmad et al., 2010). The interplay between RH and aerosol optical properties has long been recognized in marine environments (Shettle and Fenn, 1979; Gerber, 1985). However, these relationships may not apply to inland and nearshore coastal waters having a weak marine influence due to the major impact of land processes (e.g., evapotranspiration (Raoufi and Beighley, 2017)) on RH. Unlike empirical models, theoretical A2010 models use a suite of AERONET measurements collected across several small open-ocean islands and three stations along the Chesapeake Bay region to constrain the aerosol properties in marine environments suited for heritage ocean-color observations (open water and coastal, respectively). Further, A2010 models have two main characteristics: 1) they have been developed by assuming two size distribution modes (i.e., fine and coarse) having an extreme range of contributions in terms of number density, and 2) they
<figure><img src="image_2.png" /><figcaption>**FIGURE 1** | Single Scattering Albedo (SSA) and asymmetry parameter (g) for selected A2010 (RH = 80%). The models are sorted according to size fractions with increasing coarse modes from A1 through A9.</figcaption></figure>

only simulate weakly absorbing aerosols by using spherical and relatively small-sized particles (i.e., dust is not included). As an example, SSA and asymmetric properties of A2010 models corresponding to RH = 80% are shown in **Figure 1**. For NASA's standard AC (Mobley et al., 2016), the selection of these aerosol models is constrained by RH sourced from assimilated ancillary data. A recent study has reported large uncertainties in RH over continental regions (Mélin et al., 2022), likely leading to incorrect selection of aerosol models and subsequent effects on derived $R_{rs}$ over inland and nearshore coastal waters.

In general, the use of empirical aerosol models (i.e., entirely derived from *in situ* measurements) for AC has been customary in terrestrial remote sensing studies (Vermote et al., 1997; Dubovik et al., 2002a; Roger et al., 2022). These models are usually calibrated with data obtained from well-known aerosol sources (e.g., biomass burning) (Sayer et al., 2014). Although several statistical techniques have been proposed for deriving global aerosol classes from airborne and land-based (e.g., AERONET) measurements (Omar et al., 2005; Giles et al., 2012; Russell et al., 2014), none of these investigations have been applied to the study of aquatic ecosystems. In essence, no comparisons between empirical and theoretical aerosol models have been made to evaluate the advantage and disadvantages of each model for representing aerosol mixtures over aquatic environments having a weak (e.g., estuaries) or strong (e.g., lake) continental influence.

Our objective is to develop and validate globally applicable empirical aerosol mixtures that can be applied to continental waters, complement A2010 models, and improve AC methods over inland and nearshore coastal waters. Our empirical models are constructed according to AERONET measurements associated with relatively low aerosol optical depth (AOD ($0.44 \le 0.4$)), the range where AC is oftentimes carried out in practice. The methodology for this contribution is organized in four subsections 1) the classification of aerosols based on clustering of aerosol properties as derived from global AERONET measurements, 2) the validation of aerosol types derived from the proposed clustering technique (i.e., Partition Around Medoids; PAM), 3) the comparisons between PAM- and

A2010-based aerosol models for microphysical properties, and 4) the demonstration of the validity and efficacy of PAM models in the context of AC.

The first subsection involves the dimensionality reduction of aerosol parameters and describes the criteria for selecting the maximum number of clusters. The identification of aerosol sources based on literature, the definition of PAM clusters in terms of aerosol types, and the presentation of the classification performance are explained in the second subsection. The verification of realistic representations of aerosol types was centered around two techniques (the dominant clusters and the Multidimensional Euclidean Distance; MED). The third subsection highlights similarities and differences between PAM and A2010 models and the benefits and caveats of using them in continental- and marine-dominated aerosol environments. Lastly, we describe experimental machine-learning (ML) AC processors to showcase the advantages of integrating PAM aerosol models for retrieving $R_{rs}$ from OLI and validating them against matchups from selected AERONET sites equipped with ocean-color radiometers (i.e., AERONET-OC).

# DATASETS AND METHODS

## AERONET Measurements

Optical and microphysical aerosol properties for the total atmospheric column were obtained from daily Version 3, Level 2.0 inversion AERONET products between June 1993 and September 2019. The geographic distribution of measurements is shown in **Figure 2**. Our original database has 74,627 records from 806 monitoring locations. The quality control of Version 3, Level 2.0 products includes an automatic cloud screening and instrument anomaly detection (e.g., obstruction of field-of-view) (Smirnov et al., 2000; Giles et al., 2019). The inversion of aerosol properties relies on fitting radiance measurements (i.e., direct Sun and angular distribution of diffuse contributions) to radiance estimates derived from a radiative transfer model (Dubovik and King, 2000). Consistently, all Sun and sky radiance measurements made by AERONET since 1993 have been performed with
<figure><img src="image_3.png" /><figcaption>**FIGURE 1** | Single Scattering Albedo (SSA) and asymmetry parameter (g) for selected A2010 (RH = 80%). The models are sorted according to size fractions with increasing coarse modes from A1 through A9.</figcaption></figure>

only simulate weakly absorbing aerosols by using spherical and relatively small-sized particles (i.e., dust is not included). As an example, SSA and asymmetric properties of A2010 models corresponding to RH = 80% are shown in **Figure 1**. For NASA's standard AC (Mobley et al., 2016), the selection of these aerosol models is constrained by RH sourced from assimilated ancillary data. A recent study has reported large uncertainties in RH over continental regions (Mélin et al., 2022), likely leading to incorrect selection of aerosol models and subsequent effects on derived $R_{rs}$ over inland and nearshore coastal waters.

In general, the use of empirical aerosol models (i.e., entirely derived from *in situ* measurements) for AC has been customary in terrestrial remote sensing studies (Vermote et al., 1997; Dubovik et al., 2002a; Roger et al., 2022). These models are usually calibrated with data obtained from well-known aerosol sources (e.g., biomass burning) (Sayer et al., 2014). Although several statistical techniques have been proposed for deriving global aerosol classes from airborne and land-based (e.g., AERONET) measurements (Omar et al., 2005; Giles et al., 2012; Russell et al., 2014), none of these investigations have been applied to the study of aquatic ecosystems. In essence, no comparisons between empirical and theoretical aerosol models have been made to evaluate the advantage and disadvantages of each model for representing aerosol mixtures over aquatic environments having a weak (e.g., estuaries) or strong (e.g., lake) continental influence.

Our objective is to develop and validate globally applicable empirical aerosol mixtures that can be applied to continental waters, complement A2010 models, and improve AC methods over inland and nearshore coastal waters. Our empirical models are constructed according to AERONET measurements associated with relatively low aerosol optical depth (AOD ($0.44 \le 0.4$)), the range where AC is oftentimes carried out in practice. The methodology for this contribution is organized in four subsections 1) the classification of aerosols based on clustering of aerosol properties as derived from global AERONET measurements, 2) the validation of aerosol types derived from the proposed clustering technique (i.e., Partition Around Medoids; PAM), 3) the comparisons between PAM- and

A2010-based aerosol models for microphysical properties, and 4) the demonstration of the validity and efficacy of PAM models in the context of AC.

The first subsection involves the dimensionality reduction of aerosol parameters and describes the criteria for selecting the maximum number of clusters. The identification of aerosol sources based on literature, the definition of PAM clusters in terms of aerosol types, and the presentation of the classification performance are explained in the second subsection. The verification of realistic representations of aerosol types was centered around two techniques (the dominant clusters and the Multidimensional Euclidean Distance; MED). The third subsection highlights similarities and differences between PAM and A2010 models and the benefits and caveats of using them in continental- and marine-dominated aerosol environments. Lastly, we describe experimental machine-learning (ML) AC processors to showcase the advantages of integrating PAM aerosol models for retrieving $R_{rs}$ from OLI and validating them against matchups from selected AERONET sites equipped with ocean-color radiometers (i.e., AERONET-OC).

# DATASETS AND METHODS

## AERONET Measurements

Optical and microphysical aerosol properties for the total atmospheric column were obtained from daily Version 3, Level 2.0 inversion AERONET products between June 1993 and September 2019. The geographic distribution of measurements is shown in **Figure 2**. Our original database has 74,627 records from 806 monitoring locations. The quality control of Version 3, Level 2.0 products includes an automatic cloud screening and instrument anomaly detection (e.g., obstruction of field-of-view) (Smirnov et al., 2000; Giles et al., 2019). The inversion of aerosol properties relies on fitting radiance measurements (i.e., direct Sun and angular distribution of diffuse contributions) to radiance estimates derived from a radiative transfer model (Dubovik and King, 2000). Consistently, all Sun and sky radiance measurements made by AERONET since 1993 have been performed with
<figure><img src="image_3.png" /><figcaption>**FIGURE 2** | The global network of AERONET sites. The sites with > 50 daily records are labeled with numbers (see **Supplementary Appendedix Table A2** for specific site names). The green and magenta circles denote the sites whose observations were utilized in clustering. The sites without SSA/AE<sub>a</sub> inversions are marked with red circles. The data for developing PAM models were collected between October-1993 and September-2019.</figcaption></figure>

<table><caption>**TABLE 1** | List of acronyms.</caption><thead><tr><th>Abbreviation</th><th>Definition</th><th>Units</th></tr></thead><tbody><tr><td>AE<sub>a</sub></td><td>Ångtrom exponent for absorption</td><td>Dimensionless</td></tr><tr><td>AE<sub>e</sub></td><td>Ångtrom exponent for extinction</td><td>Dimensionless</td></tr><tr><td>AOD</td><td>Aerosol Optical Depth</td><td>Dimensionless</td></tr><tr><td>C<sub>n</sub></td><td>Number concentration</td><td>Particles µm<sup>-3</sup></td></tr><tr><td>CS</td><td>Relative cluster size</td><td>%</td></tr><tr><td>CSG</td><td>Cluster structure goodness</td><td>%</td></tr><tr><td>C<sub>v</sub></td><td>Volume concentration</td><td>µm<sup>3</sup> µm<sup>-2</sup></td></tr><tr><td>FC<sub>v</sub><sup>fine</sup></td><td>Fine-mode contribution to C<sub>v</sub></td><td>%</td></tr><tr><td>g</td><td>Asymmetry parameter</td><td>Dimensionless</td></tr><tr><td>K<sub>ext</sub></td><td>Extinction coefficient</td><td>m<sup>-1</sup></td></tr><tr><td>MED<sup>min</sup></td><td>Minimum multidimensional Euclidean distance</td><td>Dimensionless</td></tr><tr><td>NDF</td><td>Number density fraction</td><td>Dimensionless</td></tr><tr><td>n<sub>i</sub></td><td>Imaginary part of refractive index</td><td>Dimensionless</td></tr><tr><td>n<sub>r</sub></td><td>Real part of refractive index</td><td>Dimensionless</td></tr><tr><td>P<sub>a</sub></td><td>Phase scattering function</td><td>sr<sup>-1</sup></td></tr><tr><td>PAM</td><td>Partition Around Medoids</td><td>NA</td></tr><tr><td>PDF</td><td>Probability Distribution Function</td><td>Dimensionless</td></tr><tr><td>PSD</td><td>Particle Size Distribution</td><td>dV/dlnr or dN/dlnr</td></tr><tr><td>θs</td><td>Scattering angle</td><td>Degrees</td></tr><tr><td>r<sup>2</sup></td><td>Coefficient of determination</td><td>Dimensionless</td></tr><tr><td>RAA</td><td>Relative Azimuth Angle</td><td>Degrees</td></tr><tr><td>R<sub>rs</sub></td><td>Remote sensing reflectance</td><td>sr<sup>-1</sup></td></tr><tr><td>ρ<sub>t</sub></td><td>Top of atmosphere reflectance</td><td>Dimensionless</td></tr><tr><td>RH</td><td>Relative Humidity</td><td>%</td></tr><tr><td>r<sub>n</sub></td><td>Arithmetic average of the radius in number density space</td><td>µm</td></tr><tr><td>σ<sub>v</sub></td><td>Broadness of PSD in number density space</td><td>ln(µm)</td></tr><tr><td>SSA</td><td>Single Scattering Albedo</td><td>Dimensionless</td></tr><tr><td>σ<sub>θ</sub></td><td>Broadness of PSD in volume space</td><td>Dimensionless</td></tr><tr><td>SZA</td><td>Solar Zenith Angle</td><td>Degrees</td></tr><tr><td>VZA</td><td>View Zenith Angle</td><td>Degrees</td></tr></tbody></table>
<figure><img src="image_4.png" /><figcaption>**FIGURE 2** | The global network of AERONET sites. The sites with > 50 daily records are labeled with numbers (see **Supplementary Appendedix Table A2** for specific site names). The green and magenta circles denote the sites whose observations were utilized in clustering. The sites without SSA/AE<sub>a</sub> inversions are marked with red circles. The data for developing PAM models were collected between October-1993 and September-2019.</figcaption></figure>

<table><caption>**TABLE 1** | List of acronyms.</caption><thead><tr><th>Abbreviation</th><th>Definition</th><th>Units</th></tr></thead><tbody><tr><td>AE<sub>a</sub></td><td>Ångtrom exponent for absorption</td><td>Dimensionless</td></tr><tr><td>AE<sub>e</sub></td><td>Ångtrom exponent for extinction</td><td>Dimensionless</td></tr><tr><td>AOD</td><td>Aerosol Optical Depth</td><td>Dimensionless</td></tr><tr><td>C<sub>n</sub></td><td>Number concentration</td><td>Particles µm<sup>-3</sup></td></tr><tr><td>CS</td><td>Relative cluster size</td><td>%</td></tr><tr><td>CSG</td><td>Cluster structure goodness</td><td>%</td></tr><tr><td>C<sub>v</sub></td><td>Volume concentration</td><td>µm<sup>3</sup> µm<sup>-2</sup></td></tr><tr><td>FC<sub>v</sub><sup>fine</sup></td><td>Fine-mode contribution to C<sub>v</sub></td><td>%</td></tr><tr><td>g</td><td>Asymmetry parameter</td><td>Dimensionless</td></tr><tr><td>K<sub>ext</sub></td><td>Extinction coefficient</td><td>m<sup>-1</sup></td></tr><tr><td>MED<sup>min</sup></td><td>Minimum multidimensional Euclidean distance</td><td>Dimensionless</td></tr><tr><td>NDF</td><td>Number density fraction</td><td>Dimensionless</td></tr><tr><td>n<sub>i</sub></td><td>Imaginary part of refractive index</td><td>Dimensionless</td></tr><tr><td>n<sub>r</sub></td><td>Real part of refractive index</td><td>Dimensionless</td></tr><tr><td>P<sub>a</sub></td><td>Phase scattering function</td><td>sr<sup>-1</sup></td></tr><tr><td>PAM</td><td>Partition Around Medoids</td><td>NA</td></tr><tr><td>PDF</td><td>Probability Distribution Function</td><td>Dimensionless</td></tr><tr><td>PSD</td><td>Particle Size Distribution</td><td>dV/dlnr or dN/dlnr</td></tr><tr><td>θs</td><td>Scattering angle</td><td>Degrees</td></tr><tr><td>r<sup>2</sup></td><td>Coefficient of determination</td><td>Dimensionless</td></tr><tr><td>RAA</td><td>Relative Azimuth Angle</td><td>Degrees</td></tr><tr><td>R<sub>rs</sub></td><td>Remote sensing reflectance</td><td>sr<sup>-1</sup></td></tr><tr><td>ρ<sub>t</sub></td><td>Top of atmosphere reflectance</td><td>Dimensionless</td></tr><tr><td>RH</td><td>Relative Humidity</td><td>%</td></tr><tr><td>r<sub>n</sub></td><td>Arithmetic average of the radius in number density space</td><td>µm</td></tr><tr><td>σ<sub>v</sub></td><td>Broadness of PSD in number density space</td><td>ln(µm)</td></tr><tr><td>SSA</td><td>Single Scattering Albedo</td><td>Dimensionless</td></tr><tr><td>σ<sub>θ</sub></td><td>Broadness of PSD in volume space</td><td>Dimensionless</td></tr><tr><td>SZA</td><td>Solar Zenith Angle</td><td>Degrees</td></tr><tr><td>VZA</td><td>View Zenith Angle</td><td>Degrees</td></tr></tbody></table>
CIMEL sun-photometers (Holben et al., 1998). The spectral configuration of these instruments is based on the ultraviolet (UV) (λ = 0.34 and 0.38 µm), visible (λ = 0.44, 0.5, 0.675 µm), as well as NIR and SWIR (λ = 0.87, 1.02 and 1.64 µm) channels that had been slightly modified in some monitoring sites (Barreto et al., 2016). The daily arithmetic average of aerosol optical inversions is reported at four wavelengths (λ = 0.44, 0.675,0.87, and 1.02 µm). It comprises microphysical (e.g., real (n<sub>r</sub>) and imaginary (n<sub>i</sub>) part of refractive index, respectively) and optical (e.g., SSA) properties (Table 1).

For the spectral range of 0.44–0.87 µm and AOD(0.44) values between 0.4 and 1, the average uncertainty for AOD, SSA, the asymmetry parameter (g), n<sub>r</sub> and n<sub>i</sub>, at λ = 0.44 µm is ±0.02, ±0.03, ±4%, ±0.04, and ±30–50%, respectively (Dubovik et al., 2002b; Johnson et al., 2009; Giles et al., 2012). Uncertainties in SSA, n<sub>r</sub>, and n<sub>i</sub> values increase for AOD(0.44) < 0.4, leading to maximum uncertainties of ±0.09, ±0.05, and ±100%, respectively (Dubovik et al., 2002a; Dubovik et al., 2002b). Note that errors on SSA retrievals, as inferred from their standard deviations, vary spectrally and with respect to aerosol types (e.g., ±0.055 and ±0.11 at λ = 1.02 µm and for dust and urban/industrial aerosol sources, respectively) (Sinyuk et al., 2020). For the visible-NIR region, the maximum uncertainty associated with the extinction and absorption Ångtrom exponents (AE<sub>e</sub> and AE<sub>a</sub>, respectively) is approximately ±9% (Wagner and Silva, 2008; Giles et al., 2012). Uncertainties in microphysical properties r<sub>v</sub>, σ<sub>v</sub> and C<sub>v</sub> (i.e., the logarithmic volume median radius, the standard deviation, and the particle volume concentration, respectively) are less influenced by AOD and may vary between ±10% and ±35% for an intermediate size range (i.e., particle radius, r = 0.1–7 µm) and up to 100% near the edges of the size distributions (Dubovik et al., 2000). These uncertainties on retrieving aerosol size properties also depend on the aerosol type (Giles et al., 2012).

Selection of Clustering Variables

The number and type (i.e., optical vs. microphysical) of aerosol
properties used for discriminating different aerosol mixtures
based on AERONET measurements vary across various studies
due, in part, to differences in the number of records and their
temporal/spatial characteristics. Here, the selection of aerosol
properties was according to five criteria: 1) minimum cross-
covariation between variables, 2) aerosol properties linked to
PSD, particle shape and absorption characteristics (e.g., AE<sub>e</sub>, g,
and SSA, respectively), 3) optical parameters having spectral
information (e.g., AE<sub>e</sub>, AE<sub>a</sub>), 4) aerosol attributes with
relatively small uncertainties (i.e., < 30%), and 5) number of
parameters leading to the maximum number of clusters (K<sub>max</sub>)
(see section 2.3). Three aerosol properties met these criteria:
SSA(0.44), AE<sub>e</sub>, and AE<sub>a</sub> (Supplementary Appendix Table A1).
Indeed, preliminary results suggested minimum covariation
between the variables above, as inferred from the coefficient of
determination (r<sup>2</sup> < 0.25). For instance, SSA(0.87) was excluded
from the analysis since preliminary results indicated a lower
number of clusters when this property was included. Spurious
AE<sub>a</sub> estimates were removed as they represented < 0.5% of the
total records used for clustering. An important assumption in our
clustering is that aerosol properties excluded as input in the PAM

classification are also unique for each cluster they are associated
with because the considered aerosol properties are substantially
correlated with those excluded. The AERONET dataset was
stratified in different ranges of AOD as suggested in several
studies (Levy et al., 2007; Kahn et al., 2016). The main
application of our aerosol classification is to identify aerosol
models applicable to AC over inland and coastal waters.
Naturally, the AC and subsequent inversions of in-water
optical properties are relatively confidently conducted if the
aerosol loading were relatively low (i.e., AOD(0.44) ≤0.4).
Thus, only aerosol measurements obtained under this
atmospheric turbidity were used for the PAM classification,
reducing the number of AERONET samples by 87%, leading
to a total of 9,352 samples. This data reduction did not
compromise the number of computed clusters, nor the general
spatial coverage and seasonal variability of aerosol properties
summarized in the original dataset (i.e., with all the AOD values
included). Although accurate SSA and AEa estimates can only be
achieved for both AOD(0.44)≥ 0.4 and large zenith angles (i.e., >
50°), we adopted daily averages of relatively low AODs that
contain instances of AERONET Level 2.0 records with AODs
exceeding the threshold. By using this AERONET Level 2.0
subset, we assume that the absorption-related properties can
be applied to the lower AOD data for the entire day. This
assumption was validated by comparing changes in mean
SSA(0.44) values as a function of mean AOD(0.44) for
different aerosol loading intervals. For example, at the GSFC
site, the coefficient of determination was close to zero suggesting
no relationships for daily mean comparisons of SSA-AOD and
AEa-AOD obtained for low (≤0.4) and intermediate (i.e., ≤ 1)
daily mean AOD values (see Supplementary Appendix
Figure A1).

Note that high atmospheric turbidity (AOD(0.44) ≥ 1) represented 15% of the AERONET dataset assembled in this study. The above pre-sorted database (N = 9,352) containing SSA and AE<sub>a</sub> inversions was split into two batches (i.e., 75% training and 25% validation data chosen based on odd/even indices (Omar et al., 2005)). Lastly, these data subsets were stratified according to AOD(0.44) thresholds (i.e., only records with AOD(0.44) ≤0.4 were preserved) before applying the PAM clustering technique (see Section 2.3 and Section 2.4).

Clustering Algorithm

The classification of aerosol properties has been investigated using unsupervised and supervised algorithms with different skills and limitations (Omar et al., 2005; Giles et al., 2012; Russell et al., 2014; Hamill et al., 2016). Here, discrete aerosol classes were estimated using an unsupervised classification based on PAM (Kaufman and Rousseeuw, 1990) (Figure 3).

This top-down clustering algorithm minimizes the distance of
each point *i* with respect to the nearest centroid or “medoid” and
requires a first guess in terms of K<sub>max</sub> (the maximum number of
clusters). Unlike the k-means algorithm (MacQueen, 1967), the
“medoid” is a data point, and PAM can be used with non-
Euclidean distances and non-spherical clusters. The
Mahalanobis Distance (MD) was chosen for the PAM analysis
since it is not sensitive to outliers and weak cross-correlations
CIMEL sun-photometers (Holben et al., 1998). The spectral configuration of these instruments is based on the ultraviolet (UV) (λ = 0.34 and 0.38 µm), visible (λ = 0.44, 0.5, 0.675 µm), as well as NIR and SWIR (λ = 0.87, 1.02 and 1.64 µm) channels that had been slightly modified in some monitoring sites (Barreto et al., 2016). The daily arithmetic average of aerosol optical inversions is reported at four wavelengths (λ = 0.44, 0.675,0.87, and 1.02 µm). It comprises microphysical (e.g., real (n<sub>r</sub>) and imaginary (n<sub>i</sub>) part of refractive index, respectively) and optical (e.g., SSA) properties (Table 1).

For the spectral range of 0.44–0.87 µm and AOD(0.44) values between 0.4 and 1, the average uncertainty for AOD, SSA, the asymmetry parameter (g), n<sub>r</sub> and n<sub>i</sub>, at λ = 0.44 µm is ±0.02, ±0.03, ±4%, ±0.04, and ±30–50%, respectively (Dubovik et al., 2002b; Johnson et al., 2009; Giles et al., 2012). Uncertainties in SSA, n<sub>r</sub>, and n<sub>i</sub> values increase for AOD(0.44) < 0.4, leading to maximum uncertainties of ±0.09, ±0.05, and ±100%, respectively (Dubovik et al., 2002a; Dubovik et al., 2002b). Note that errors on SSA retrievals, as inferred from their standard deviations, vary spectrally and with respect to aerosol types (e.g., ±0.055 and ±0.11 at λ = 1.02 µm and for dust and urban/industrial aerosol sources, respectively) (Sinyuk et al., 2020). For the visible-NIR region, the maximum uncertainty associated with the extinction and absorption Ångtrom exponents (AE<sub>e</sub> and AE<sub>a</sub>, respectively) is approximately ±9% (Wagner and Silva, 2008; Giles et al., 2012). Uncertainties in microphysical properties r<sub>v</sub>, σ<sub>v</sub> and C<sub>v</sub> (i.e., the logarithmic volume median radius, the standard deviation, and the particle volume concentration, respectively) are less influenced by AOD and may vary between ±10% and ±35% for an intermediate size range (i.e., particle radius, r = 0.1–7 µm) and up to 100% near the edges of the size distributions (Dubovik et al., 2000). These uncertainties on retrieving aerosol size properties also depend on the aerosol type (Giles et al., 2012).

## Selection of Clustering Variables

The number and type (i.e., optical vs. microphysical) of aerosol properties used for discriminating different aerosol mixtures based on AERONET measurements vary across various studies due, in part, to differences in the number of records and their temporal/spatial characteristics. Here, the selection of aerosol properties was according to five criteria: 1) minimum cross-covariation between variables, 2) aerosol properties linked to PSD, particle shape and absorption characteristics (e.g., AE<sub>e</sub>, g, and SSA, respectively), 3) optical parameters having spectral information (e.g., AE<sub>e</sub>, AE<sub>a</sub>), 4) aerosol attributes with relatively small uncertainties (i.e., < 30%), and 5) number of parameters leading to the maximum number of clusters (K<sub>max</sub>) (see section 2.3). Three aerosol properties met these criteria: SSA(0.44), AE<sub>e</sub>, and AE<sub>a</sub> (Supplementary Appendix Table A1). Indeed, preliminary results suggested minimum covariation between the variables above, as inferred from the coefficient of determination (r<sup>2</sup> < 0.25). For instance, SSA(0.87) was excluded from the analysis since preliminary results indicated a lower number of clusters when this property was included. Spurious AE<sub>a</sub> estimates were removed as they represented < 0.5% of the total records used for clustering. An important assumption in our clustering is that aerosol properties excluded as input in the PAM

classification are also unique for each cluster they are associated with because the considered aerosol properties are substantially correlated with those excluded. The AERONET dataset was stratified in different ranges of AOD as suggested in several studies (Levy et al., 2007; Kahn et al., 2016). The main application of our aerosol classification is to identify aerosol models applicable to AC over inland and coastal waters. Naturally, the AC and subsequent inversions of in-water optical properties are relatively confidently conducted if the aerosol loading were relatively low (i.e., AOD(0.44) ≤0.4). Thus, only aerosol measurements obtained under this atmospheric turbidity were used for the PAM classification, reducing the number of AERONET samples by 87%, leading to a total of 9,352 samples. This data reduction did not compromise the number of computed clusters, nor the general spatial coverage and seasonal variability of aerosol properties summarized in the original dataset (i.e., with all the AOD values included). Although accurate SSA and AEa estimates can only be achieved for both AOD(0.44)≥ 0.4 and large zenith angles (i.e., > 50°), we adopted daily averages of relatively low AODs that contain instances of AERONET Level 2.0 records with AODs exceeding the threshold. By using this AERONET Level 2.0 subset, we assume that the absorption-related properties can be applied to the lower AOD data for the entire day. This assumption was validated by comparing changes in mean SSA(0.44) values as a function of mean AOD(0.44) for different aerosol loading intervals. For example, at the GSFC site, the coefficient of determination was close to zero suggesting no relationships for daily mean comparisons of SSA-AOD and AEa-AOD obtained for low (≤0.4) and intermediate (i.e., ≤ 1) daily mean AOD values (see Supplementary Appendix Figure A1).

Note that high atmospheric turbidity (AOD(0.44) ≥ 1) represented 15% of the AERONET dataset assembled in this study. The above pre-sorted database (N = 9,352) containing SSA and AEa inversions was split into two batches (i.e., 75% training and 25% validation data chosen based on odd/even indices (Omar et al., 2005)). Lastly, these data subsets were stratified according to AOD(0.44) thresholds (i.e., only records with AOD(0.44) ≤0.4 were preserved) before applying the PAM clustering technique (see Section 2.3 and Section 2.4).

## Clustering Algorithm

The classification of aerosol properties has been investigated using unsupervised and supervised algorithms with different skills and limitations (Omar et al., 2005; Giles et al., 2012; Russell et al., 2014; Hamill et al., 2016). Here, discrete aerosol classes were estimated using an unsupervised classification based on PAM (Kaufman and Rousseeuw, 1990) (Figure 3).

This top-down clustering algorithm minimizes the distance of each point *i* with respect to the nearest centroid or “medoid” and requires a first guess in terms of K<sub>max</sub> (the maximum number of clusters). Unlike the k-means algorithm (MacQueen, 1967), the “medoid” is a data point, and PAM can be used with non-Euclidean distances and non-spherical clusters. The Mahalanobis Distance (MD) was chosen for the PAM analysis since it is not sensitive to outliers and weak cross-correlations
<figure><img src="image_4.png" /><figcaption>FIGURE 3 | PAM classification of aerosols.</figcaption></figure>

between variables. MD measures the dissimilarity between two vectors (u and v) with covariance S (Varmuza and Filzmoser, 2016):

$$MD(u, v) = \sqrt{(u - v)^T S^{-1} (u - v)} \quad (1)$$

where superscript T denotes the matrix transpose, and u, v vectors are composed of pre-selected aerosol properties. Kmax was estimated using the Silhouette index (Sil) (Rousseeuw, 1987) after initializing PAM with 30 potential clusters (i.e., initial guess). This initial guess corresponds to an intermediate number of aerosol mixtures with respect to the range reported in the literature (i.e., 4 to 74 aerosol models) (Omar et al., 2005; Giles et al., 2012; Russell et al., 2014; Li et al., 2015; Hamill et al., 2016)). Sil is defined as:

$$Sil^k(m) = \frac{(b(m) - a(m))}{max\{a(m), b(m)\}} \quad (2)$$

where b(m) and a(m) are the minimum dissimilarities ("distance") from sample m to samples not belonging to the cluster k. The Silhouette magnitude varies between -1 and 1 and measures how similar an object is relative to its own cluster with values near +1 (cohesion) compared to other groups with values near -1 (separation). The criterion here for estimating the final Kmax used in each clustering process was a Sil greater or equal to 0.51 (i.e., a threshold above which the clustering structure is at least "reasonable") (Kaufman and Rousseeuw, 1990). The remaining Sil-based clustering structure definitions are "no substantial", "weak", and "strong" (i.e., -1 to 0.25, 0.26-0.5, and >71, respectively). Preliminary results suggest insignificance of two small-sized clusters (N < 30) that were excluded from the analysis due to extremely low (high) SSA (AEa) values likely related to the relatively high uncertainties of these

parameters. Thus, a group of records was considered a cluster only if that potential cluster has at least 50 daily records with Sil ≥ 0.51.

The relative cluster size (CS) was defined as the number of records belonging to cluster k divided by the number of records of all clusters encompassing 9,352 observations. Lastly, the cluster structure goodness (CSG) was defined as the number of records of cluster k with Sil ≥ 0.51 divided by the total number of records in cluster k. CSG is independent of the number of measurements and indicative of how "compact" (i.e., homogeneous in terms of aerosol properties) each cluster is. No CSG thresholds were used for defining the maximum number of clusters as estimated from PAM.

## Evaluation of Clusters

Ground-truthing of aerosol models derived from clustering of AERONET measurements is commonly based on known aerosol sources (e.g., areas with continuous industrial activity) (Omar et al., 2005) or in situ airborne experiments (e.g., Aerosol Characterization Experiment-Asia) (Kim et al., 2007). An aerosol source may only be dominant within a month or season in some areas rather than during an entire year (e.g., dust in spring at Kanpur, India) (Giles et al., 2012). In general, aerosol sources are linked to five primitive aerosol types: 1. dust, 2. mixed, 3. urban/industrial, 4. biomass burning, and 5. marine, hereafter DU, MIX, U/I, BB, and MAR, respectively (Omar et al., 2005).

The first step in evaluating the clustering results was the association of each common AERONET site (i.e., used by PAM and cited in the literature) with one primitive aerosol type that was dominant throughout the year (Omar et al., 2005; Giles et al., 2012; Russell et al., 2014) (Supplementary Appendix Table A2). The second step was to label each cluster in
<figure><img src="image_5.png" /><figcaption>FIGURE 3 | PAM classification of aerosols.</figcaption></figure>

between variables. MD measures the dissimilarity between two vectors (u and v) with covariance S (Varmuza and Filzmoser, 2016):

$$MD(u, v) = \sqrt{(u - v)^T S^{-1} (u - v)} \quad (1)$$

where superscript T denotes the matrix transpose, and u, v vectors are composed of pre-selected aerosol properties. Kmax was estimated using the Silhouette index (Sil) (Rousseeuw, 1987) after initializing PAM with 30 potential clusters (i.e., initial guess). This initial guess corresponds to an intermediate number of aerosol mixtures with respect to the range reported in the literature (i.e., 4 to 74 aerosol models) (Omar et al., 2005; Giles et al., 2012; Russell et al., 2014; Li et al., 2015; Hamill et al., 2016)). Sil is defined as:

$$Sil^k(m) = \frac{(b(m) - a(m))}{max\{a(m), b(m)\}} \quad (2)$$

where b(m) and a(m) are the minimum dissimilarities (“distance”) from sample m to samples not belonging to the cluster k. The Silhouette magnitude varies between -1 and 1 and measures how similar an object is relative to its own cluster with values near +1 (cohesion) compared to other groups with values near -1 (separation). The criterion here for estimating the final Kmax used in each clustering process was a Sil greater or equal to 0.51 (i.e., a threshold above which the clustering structure is at least “reasonable”) (Kaufman and Rousseeuw, 1990). The remaining Sil-based clustering structure definitions are “no substantial”, “weak”, and “strong” (i.e., -1 to 0.25, 0.26-0.5, and >71, respectively). Preliminary results suggest insignificance of two small-sized clusters (N < 30) that were excluded from the analysis due to extremely low (high) SSA (AEa) values likely related to the relatively high uncertainties of these

parameters. Thus, a group of records was considered a cluster only if that potential cluster has at least 50 daily records with Sil ≥ 0.51.

The relative cluster size (CS) was defined as the number of records belonging to cluster k divided by the number of records of all clusters encompassing 9,352 observations. Lastly, the cluster structure goodness (CSG) was defined as the number of records of cluster k with Sil ≥ 0.51 divided by the total number of records in cluster k. CSG is independent of the number of measurements and indicative of how “compact” (i.e., homogeneous in terms of aerosol properties) each cluster is. No CSG thresholds were used for defining the maximum number of clusters as estimated from PAM.

## Evaluation of Clusters

Ground-truthing of aerosol models derived from clustering of AERONET measurements is commonly based on known aerosol sources (e.g., areas with continuous industrial activity) (Omar et al., 2005) or in situ airborne experiments (e.g., Aerosol Characterization Experiment-Asia) (Kim et al., 2007). An aerosol source may only be dominant within a month or season in some areas rather than during an entire year (e.g., dust in spring at Kanpur, India) (Giles et al., 2012). In general, aerosol sources are linked to five primitive aerosol types: 1. dust, 2. mixed, 3. urban/industrial, 4. biomass burning, and 5. marine, hereafter DU, MIX, U/I, BB, and MAR, respectively (Omar et al., 2005).

The first step in evaluating the clustering results was the association of each common AERONET site (i.e., used by PAM and cited in the literature) with one primitive aerosol type that was dominant throughout the year (Omar et al., 2005; Giles et al., 2012; Russell et al., 2014) (Supplementary Appendix Table A2). The second step was to label each cluster in
terms of “most likely” aerosol type by analyzing the frequency and presence of sites described in step 1 in each cluster. Thus, more than one cluster can be aggregated within the same aerosol type resulting in “aerosol sub-types”. The third step was the prediction of the aerosol type in each reference site as estimated from PAM and comparison with the initial label assigned in step 1. To accomplish this last task, the PAM results were validated with respect to different geographic locations and periods of the year (annual and seasonal analysis, respectively).

The annual analysis encompassed the evaluation of two metrics: 1) dominant annual clusters, and 2) minimum Multidimensional Euclidean Distance (MED<sup>min</sup>). Comparisons for MAR aerosols were not performed due to the lack of data given the relatively low AOD of marine environments and known difficulties for obtaining accurate SSA and AE<sub>a</sub> inversions. A dominant cluster and associated aerosol type at the AERONET site j were defined as the most frequent cluster if the minimum number of records was greater than two and no co-dominance was present (i.e., at least two aerosol types have the same frequency). The MED<sup>min</sup> metrics were applied to an independent AERONET dataset created in section 2.2 computed at the monthly scale based on dominant clusters derived from PAM. The MED<sup>min</sup> calculations were focused on four AERONET sites (Solar Village, XiangHe, Mexico City, and Mongu) where a consensus among studies regarding the annual dominance or co-dominance of aerosol sources (DU, MIX, U/I, and BB, respectively) were found (Omar et al., 2005; Giles et al., 2012; Russell et al., 2014; Hamill et al., 2016) (Supplementary Appendix Table A2). The MED is equivalent to the norm of the difference between u and v vectors:

$$MED^{\min} = \min \left\{ \sqrt{(u_{1,m} - v_{1,m})^2 + (u_{2,m} - v_{2,m})^2 + (u_{3,m} - v_{3,m})^2} \right\} \quad (3)$$

where u and v represent mean aerosol properties for dominant PAM-derived cluster and independent dataset, respectively, for month m, subscripts 1, 2, and 3 indicate the arithmetic average of SSA(0.44), AE<sub>e</sub>, and AE<sub>a</sub>, respectively. This arithmetic average is multi-annual for PAM-derived classes and monthly for the validation dataset. The aerosol type associated with the PAM-derived cluster with the minimum MED was selected for each site and month. The classification performance of PAM for aerosol types was quantified as unity (Number of mismatches/(Number of matches + mismatches)). For the same aerosol type, matches for dominant annual clusters were computed based on several previous study sites. Conversely, matches for MED<sup>min</sup> were examined with respect to a single reference location but using multiple monthly estimates. Also, reference aerosol sources for the MED<sup>min</sup> technique were based on PAM-derived types.

Seasonal validations were examined at three AERONET sites (i.e., GSFC, Kanpur and Banizoumbou, Figure 3, Supplementary Appendix Table A1) where the aerosol composition during different periods of the year has been investigated in detail (Omar et al., 2005; Levy et al., 2007; Giles et al., 2012; Russell et al., 2014). The analysis consisted of characterizing changes in

PAM-derived aerosol types among seasons and comparing the patterns described in published studies.

## Empirical vs. Theoretical Models

The optical inversion implemented by AERONET for retrieving aerosol optical properties differs from that applied in A2010. The AERONET radiative transfer model has additional assumptions regarding surface reflectance and multiple scattering effects (García et al., 2008; Giles et al., 2012; Sinyuk et al., 2020). Thus, the optical properties of PAM-derived and A2010 aerosol models may not be compared unless Mie simulations are performed for calculating the optical properties of PAM and A2010 models. The simulations for “adjusted” A2010 models should be done using mean volume density concentration fractions derived from each PAM cluster. To simplify the calculations, a first attempt to investigate the similarities and differences between PAM and RH-based A2010 models was made using the refractive index (n<sub>r</sub> and n<sub>i</sub>) and two microphysical properties defined in the number density space (the mean radius and width of PSD, r<sub>n</sub> and σ<sub>n</sub>, respectively, Table 1). Since AERONET inversions are carried out with PSDs in volume space, a transformation was performed to compute PSD-related parameters in number density units as in A2010. Further, unlike A2010, the refractive index in PAM-derived models is not resolved as a function of PSD modes. Thus, n<sub>r</sub> and n<sub>i</sub> values for fine and coarse modes of A2010 models were merged by using the number density function (NDF) of every RH case study (i.e., 30, 50, 70, 75, 80, 85, 90, and 95%) as a weighting factor:

$$\langle n_x \rangle = n_x^{\text{fine}} NDF^{\text{fine}} + n_x^{\text{coarse}} NDF^{\text{coarse}} \quad (4)$$

where <> indicates the NDF-weighted parameter, and x represents n<sub>r</sub> or n<sub>i</sub>. As a result of the above formulation, the PAM and RH-based aerosol mixtures can be directly compared and assessed for their complementary nature and/or differences.

## Atmospheric Correction Case Study

For a practical assessment of the efficacy of the PAM models in the context of AC, we integrated our models into an experimental AC framework. Adapting the state-of-the-art processors (Pahlevan et al., 2021a) was unachievable due to a number of reasons, including, their theoretical limitations for handling absorbing aerosols (e.g., they cannot be detected using the NIR band ratio (Chomko and Gordon, 1998)) and the difficulty in their formulation and implementation (e.g., re-generating Rayleigh look-up-tables, single-scattering to multiple-scattering conversion). Inevitably, we developed two machine learning (ML) processors to convert top-of-atmosphere reflectance (ρ<sub>t</sub>) to R<sub>rs</sub> for Landsat-8/OLI data and then evaluated the quality of retrievals over selected AERONET-OC sites. To train the ML models, a coupled ocean-atmosphere, vector radiative transfer (RT) code adapted to ingest in situ (OLI-like) R<sub>rs</sub> was employed (Zhai et al., 2009). For this demonstration, OLI imagery was considered because of its utility across global inland and nearshore coastal waters and its limited viewing zenith angles (VZA < 7.5°), reducing the required number of RT simulations.
terms of “most likely” aerosol type by analyzing the frequency and presence of sites described in step 1 in each cluster. Thus, more than one cluster can be aggregated within the same aerosol type resulting in “aerosol sub-types”. The third step was the prediction of the aerosol type in each reference site as estimated from PAM and comparison with the initial label assigned in step 1. To accomplish this last task, the PAM results were validated with respect to different geographic locations and periods of the year (annual and seasonal analysis, respectively).

The annual analysis encompassed the evaluation of two metrics: 1) dominant annual clusters, and 2) minimum Multidimensional Euclidean Distance (MED<sup>min</sup>). Comparisons for MAR aerosols were not performed due to the lack of data given the relatively low AOD of marine environments and known difficulties for obtaining accurate SSA and AE<sub>a</sub> inversions. A dominant cluster and associated aerosol type at the AERONET site j were defined as the most frequent cluster if the minimum number of records was greater than two and no co-dominance was present (i.e., at least two aerosol types have the same frequency). The MED<sup>min</sup> metrics were applied to an independent AERONET dataset created in section 2.2 computed at the monthly scale based on dominant clusters derived from PAM. The MED<sup>min</sup> calculations were focused on four AERONET sites (Solar Village, XiangHe, Mexico City, and Mongu) where a consensus among studies regarding the annual dominance or co-dominance of aerosol sources (DU, MIX, U/I, and BB, respectively) were found (Omar et al., 2005; Giles et al., 2012; Russell et al., 2014; Hamill et al., 2016) (Supplementary Appendix Table A2). The MED is equivalent to the norm of the difference between u and v vectors:

$$MED^{\min} = \min \left\{ \sqrt{(u_{1,m} - v_{1,m})^2 + (u_{2,m} - v_{2,m})^2 + (u_{3,m} - v_{3,m})^2} \right\} \quad (3)$$

where u and v represent mean aerosol properties for dominant PAM-derived cluster and independent dataset, respectively, for month m, subscripts 1, 2, and 3 indicate the arithmetic average of SSA(0.44), AE<sub>e</sub>, and AE<sub>a</sub>, respectively. This arithmetic average is multi-annual for PAM-derived classes and monthly for the validation dataset. The aerosol type associated with the PAM-derived cluster with the minimum MED was selected for each site and month. The classification performance of PAM for aerosol types was quantified as unity (Number of mismatches/(Number of matches + mismatches)). For the same aerosol type, matches for dominant annual clusters were computed based on several previous study sites. Conversely, matches for MED<sup>min</sup> were examined with respect to a single reference location but using multiple monthly estimates. Also, reference aerosol sources for the MED<sup>min</sup> technique were based on PAM-derived types.

Seasonal validations were examined at three AERONET sites (i.e., GSFC, Kanpur and Banizoumbou, Figure 3, Supplementary Appendix Table A1) where the aerosol composition during different periods of the year has been investigated in detail (Omar et al., 2005; Levy et al., 2007; Giles et al., 2012; Russell et al., 2014). The analysis consisted of characterizing changes in

PAM-derived aerosol types among seasons and comparing the patterns described in published studies.

## Empirical vs. Theoretical Models

The optical inversion implemented by AERONET for retrieving aerosol optical properties differs from that applied in A2010. The AERONET radiative transfer model has additional assumptions regarding surface reflectance and multiple scattering effects (García et al., 2008; Giles et al., 2012; Sinyuk et al., 2020). Thus, the optical properties of PAM-derived and A2010 aerosol models may not be compared unless Mie simulations are performed for calculating the optical properties of PAM and A2010 models. The simulations for “adjusted” A2010 models should be done using mean volume density concentration fractions derived from each PAM cluster. To simplify the calculations, a first attempt to investigate the similarities and differences between PAM and RH-based A2010 models was made using the refractive index (n<sub>r</sub> and n<sub>i</sub>) and two microphysical properties defined in the number density space (the mean radius and width of PSD, r<sub>n</sub> and σ<sub>n</sub>, respectively, Table 1). Since AERONET inversions are carried out with PSDs in volume space, a transformation was performed to compute PSD-related parameters in number density units as in A2010. Further, unlike A2010, the refractive index in PAM-derived models is not resolved as a function of PSD modes. Thus, n<sub>r</sub> and n<sub>i</sub> values for fine and coarse modes of A2010 models were merged by using the number density function (NDF) of every RH case study (i.e., 30, 50, 70, 75, 80, 85, 90, and 95%) as a weighting factor:

$$\langle n_x \rangle = n_x^{\text{fine}} NDF^{\text{fine}} + n_x^{\text{coarse}} NDF^{\text{coarse}} \quad (4)$$

where <> indicates the NDF-weighted parameter, and x represents n<sub>r</sub> or n<sub>i</sub>. As a result of the above formulation, the PAM and RH-based aerosol mixtures can be directly compared and assessed for their complementary nature and/or differences.

## Atmospheric Correction Case Study

For a practical assessment of the efficacy of the PAM models in the context of AC, we integrated our models into an experimental AC framework. Adapting the state-of-the-art processors (Pahlevan et al., 2021a) was unachievable due to a number of reasons, including, their theoretical limitations for handling absorbing aerosols (e.g., they cannot be detected using the NIR band ratio (Chomko and Gordon, 1998)) and the difficulty in their formulation and implementation (e.g., re-generating Rayleigh look-up-tables, single-scattering to multiple-scattering conversion). Inevitably, we developed two machine learning (ML) processors to convert top-of-atmosphere reflectance (ρ<sub>t</sub>) to R<sub>rs</sub> for Landsat-8/OLI data and then evaluated the quality of retrievals over selected AERONET-OC sites. To train the ML models, a coupled ocean-atmosphere, vector radiative transfer (RT) code adapted to ingest in situ (OLI-like) R<sub>rs</sub> was employed (Zhai et al., 2009). For this demonstration, OLI imagery was considered because of its utility across global inland and nearshore coastal waters and its limited viewing zenith angles (VZA < 7.5°), reducing the required number of RT simulations.
The input to the RT code consisted of in situ $R_{rs}$ (N = 215) together with a reasonable range of solar zenith angles (SZA) (20–60° with increments of 10°) and aerosol optical depths at 550 nm (0.03–0.23 with steps of 0.04). The two ML processors differed only in terms of the aerosol models used. One processor (hereafter ACI) applied simulation results (paired sets of $R_{rs}$ and $\rho_t$) generated by incorporating A2010 aerosol models into the RT code. The other processor (hereafter ACII) was developed using RT simulations supplied with both A2010 and PAM models. Comparing the performances of ACI and ACII elucidates the effectiveness of PAM models when used in conjunction with A2010 models. Table 5 provides the information on the aerosol size distribution (effective radius and variance) and refractive index at 0.44 and 0.87 µm for each of the 27 aerosol models, which can be used to obtain the refractive index at an arbitrary wavelength through linear interpolation. This information is used by the Lorenz–Mie theory to calculate SSA and phase matrix, which is then fed to the radiative transfer model for multiple scattering calculations. AOD is a product of aerosol column number density and extinction cross-section, which is used as an independent input to the RT model. The A2010 models were incorporated into the model in a similar fashion.

To further constrain the number of simulations, only a subset of A2010 models representing fine-mode, weakly-to-non absorbing models (indexed 1 through 7) modulated by a few RH cases (50, 70, and 80%) were applied (we assumed an absence of oceanic aerosols). This choice of A2010 models was determined from Landsat-8/OLI matchups processed via the SeaWiFS Data Analysis System (SeaDAS). The final combination of aerosols for ACI yielded 21 (= 7 × 3) different aerosol models and 48 (=27 + 21) for ACII. The rest of the input parameters, including atmospheric gaseous concentrations (NO₂ and O₃ at 0.17 and 290 Dobson, respectively), water vapor amount (3 g cm⁻²), wind speed (3 m s⁻¹), and surface pressure (1,017 mb), were assumed constant.

Our ML model is a class of neural networks referred to as Mixture Density Networks (MDNs) (Bishop, 1994), which learn conditional probabilities for input ($\rho_t$) and output ($R_{rs}$) variables and are capable of modeling multimodal distributions in the target space. Due to these characteristics, MDNs proved adept at tackling the non-unique nature of the inverse problem common to remote sensing science and applications, including the aquatic remote sensing field (Smith et al., 2021). This inherent trait separates our experimental AC processors from previously designed ML models for AC (Brajard et al., 2006; Schroeder et al., 2007; Fan et al., 2021). Valid Level 2 AERONET-OC radiometric products ($R_{rs}$) matching cloud-free OLI observations [selected according to SeaDAS exclusion criteria (Werdell and Bailey, 2005)] were utilized for assessing the performance of ACI and ACII. An adequate number of matchup datasets (N = 173) was extracted mainly from coastal sites, such as LISCO (N = 7), Lake Erie (N = 8), MVCO (N = 24), Venise (N = 94), and Wave CIS (N = 40) representative of blue and/or green waters. Following a preliminary evaluation of the shapes and magnitudes of these AERONET-OC $R_{rs}$ products, a small subset (N = 215) of an in situ $R_{rs}$ database (Balasubramanian et al., 2020; Pahlevan et al., 2020; Smith

et al., 2021) resampled with OLI's relative spectral response functions (Supplementary Appendix Figure A2) was identified and incorporated as input to the RT code. It is recalled that the development of these two processors is merely intended for this demonstration, and they are not expected to operate across global aquatic ecosystems.

More specific information that corresponds to the development of MDNs is enlisted in Table 2. The input and output features were log-transformed and scaled with mean and standard deviations. Like previous studies (Smith et al., 2021; O'Shea et al., 2021), the final model prediction is the median of 10 ensembles (i.e., MDN iterations). For the matchup assessments, median $\rho_t$ computed from 3 x 3-element windows centered around AERONET-OC sites were chosen to represent OLI observations. OLI $\rho_t$ spectra were further normalized by the Earth-Sun distance squared to eliminate the associated dependency on this variability. Clouds, cloud shadows, and hazy matchups were identified and excluded using SeaDAS criteria (e.g., $\rho_t$ (2201) > 0.018) (Pahlevan et al., 2017b). Note that, in this exercise, no spectral band adjustments were carried out, due to the comparative nature of the analysis (ACI vs. ACII).

For performance evaluations, we applied the metrics below

$$X = \log_{10}(\hat{Q}/Q)$$

$$MdSA = 100 \times (10^{Y1} - 1) [\%] \text{where } Y1 = \text{Median} |X| \quad (5)$$

$$MSA = 100 \times (10^{Y2} - 1) [\%] \text{where } Y2 = \text{Mean} |X| \quad (6)$$

$$RMSLE = \sqrt{\frac{\sum_{1}^{N} X^2}{N}} \quad (7)$$

where $\hat{Q}$ and Q correspond to satellite-derived $R_{rs}$ and in situ $R_{rs}$, respectively, N is the number of AERONET-OC matchups, MdSA stands for median symmetric accuracy, MSA is mean symmetric accuracy, and RMLE represents root-mean-squared logarithmic error (Morley et al., 2018). These metrics are robust and capture different aspects of model performances (Pahlevan et al., 2021a).

# RESULTS

## Clusters and Associated Aerosol Properties

The unsupervised classification of aerosols based on PAM resulted in 27 clusters (Table 3). This number of aerosol mixtures is relative as it may vary for higher or lower Silhouette thresholds (e.g., 1 and 200 mixtures, respectively). However, an intermediate Sil index (≥0.51) was adopted to reduce the processing time without compromising the number of potential aerosol classes reported in the literature (Omar et al., 2005; Giles et al., 2012; Russell et al., 2014; Hamill et al., 2016). Each PAM cluster represents an aerosol model mixture with specific microphysical and radiative properties. Also, each cluster is composed of multiple AERONET sites and measurements obtained during different periods of the year. To interpret regional and temporal changes in the aerosol compositions, the PAM models were assigned to five broad aerosol types based on AERONET sites widely cited as aerosol sources
The input to the RT code consisted of *in situ* $R_{rs}$ (N = 215) together with a reasonable range of solar zenith angles (SZA) (20–60° with increments of 10°) and aerosol optical depths at 550 nm (0.03–0.23 with steps of 0.04). The two ML processors differed only in terms of the aerosol models used. One processor (hereafter ACI) applied simulation results (paired sets of $R_{rs}$ and $\rho_t$) generated by incorporating A2010 aerosol models into the RT code. The other processor (hereafter ACII) was developed using RT simulations supplied with both A2010 and PAM models. Comparing the performances of ACI and ACII elucidates the effectiveness of PAM models when used in conjunction with A2010 models. Table 5 provides the information on the aerosol size distribution (effective radius and variance) and refractive index at 0.44 and 0.87 µm for each of the 27 aerosol models, which can be used to obtain the refractive index at an arbitrary wavelength through linear interpolation. This information is used by the Lorenz–Mie theory to calculate SSA and phase matrix, which is then fed to the radiative transfer model for multiple scattering calculations. AOD is a product of aerosol column number density and extinction cross-section, which is used as an independent input to the RT model. The A2010 models were incorporated into the model in a similar fashion.

To further constrain the number of simulations, only a subset of A2010 models representing fine-mode, weakly-to-non absorbing models (indexed 1 through 7) modulated by a few RH cases (50, 70, and 80%) were applied (we assumed an absence of oceanic aerosols). This choice of A2010 models was determined from Landsat-8/OLI matchups processed via the SeaWiFS Data Analysis System (SeaDAS). The final combination of aerosols for ACI yielded 21 (= 7 × 3) different aerosol models and 48 (=27 + 21) for ACII. The rest of the input parameters, including atmospheric gaseous concentrations (NO₂ and O₃ at 0.17 and 290 Dobson, respectively), water vapor amount (3 g cm⁻²), wind speed (3 m s⁻¹), and surface pressure (1,017 mb), were assumed constant.

Our ML model is a class of neural networks referred to as Mixture Density Networks (MDNs) (Bishop, 1994), which learn conditional probabilities for input ($\rho_t$) and output ($R_{rs}$) variables and are capable of modeling multimodal distributions in the target space. Due to these characteristics, MDNs proved adept at tackling the non-unique nature of the inverse problem common to remote sensing science and applications, including the aquatic remote sensing field (Smith et al., 2021). This inherent trait separates our experimental AC processors from previously designed ML models for AC (Brajard et al., 2006; Schroeder et al., 2007; Fan et al., 2021). Valid Level 2 AERONET-OC radiometric products ($R_{rs}$) matching cloud-free OLI observations [selected according to SeaDAS exclusion criteria (Werdell and Bailey, 2005)] were utilized for assessing the performance of ACI and ACII. An adequate number of matchup datasets (N = 173) was extracted mainly from coastal sites, such as LISCO (N = 7), Lake Erie (N = 8), MVCO (N = 24), Venise (N = 94), and Wave CIS (N = 40) representative of blue and/or green waters. Following a preliminary evaluation of the shapes and magnitudes of these AERONET-OC $R_{rs}$ products, a small subset (N = 215) of an *in situ* $R_{rs}$ database (Balasubramanian et al., 2020; Pahlevan et al., 2020; Smith

et al., 2021) resampled with OLI's relative spectral response functions (Supplementary Appendix Figure A2) was identified and incorporated as input to the RT code. It is recalled that the development of these two processors is merely intended for this demonstration, and they are not expected to operate across global aquatic ecosystems.

More specific information that corresponds to the development of MDNs is enlisted in Table 2. The input and output features were log-transformed and scaled with mean and standard deviations. Like previous studies (Smith et al., 2021; O’Shea et al., 2021), the final model prediction is the median of 10 ensembles (i.e., MDN iterations). For the matchup assessments, median $\rho_t$ computed from 3 x 3-element windows centered around AERONET-OC sites were chosen to represent OLI observations. OLI $\rho_t$ spectra were further normalized by the Earth-Sun distance squared to eliminate the associated dependency on this variability. Clouds, cloud shadows, and hazy matchups were identified and excluded using SeaDAS criteria (e.g., $\rho_t$ (2201) > 0.018) (Pahlevan et al., 2017b). Note that, in this exercise, no spectral band adjustments were carried out, due to the comparative nature of the analysis (ACI vs. ACII).

For performance evaluations, we applied the metrics below

$$X = \log_{10}(\hat{Q}/Q) \qquad (5)$$

$$MdSA = 100 \times (10^{Y1} - 1) [\%] \text{where } Y1 = \text{Median} |X| \qquad (5)$$

$$MSA = 100 \times (10^{Y2} - 1) [\%] \text{where } Y2 = \text{Mean}|X| \qquad (6)$$

$$RMSLE = \sqrt{\frac{\sum_{1}^{N} X^2}{N}} \qquad (7)$$

where $\hat{Q}$ and Q correspond to satellite-derived $R_{rs}$ and *in situ* $R_{rs}$, respectively, N is the number of AERONET-OC matchups, MdSA stands for median symmetric accuracy, MSA is mean symmetric accuracy, and RMLE represents root-mean-squared logarithmic error (Morley et al., 2018). These metrics are robust and capture different aspects of model performances (Pahlevan et al., 2021a).

# RESULTS

## Clusters and Associated Aerosol Properties

The unsupervised classification of aerosols based on PAM resulted in 27 clusters (Table 3). This number of aerosol mixtures is relative as it may vary for higher or lower Silhouette thresholds (e.g., 1 and 200 mixtures, respectively). However, an intermediate Sil index (≥0.51) was adopted to reduce the processing time without compromising the number of potential aerosol classes reported in the literature (Omar et al., 2005; Giles et al., 2012; Russell et al., 2014; Hamill et al., 2016). Each PAM cluster represents an aerosol model mixture with specific microphysical and radiative properties. Also, each cluster is composed of multiple AERONET sites and measurements obtained during different periods of the year. To interpret regional and temporal changes in the aerosol compositions, the PAM models were assigned to five broad aerosol types based on AERONET sites widely cited as aerosol sources
<table><caption>**TABLE 2** | Properties of MDNs used for developing our experimental AC processors.</caption><thead><tr><th colspan="2">Network hyperparameters</th></tr></thead><tbody><tr><td>Training Iterations</td><td>5,000</td></tr><tr><td>Number of Gaussians</td><td>5</td></tr><tr><td>Neurons per hidden layer</td><td>30</td></tr><tr><td>Hidden layers</td><td>5</td></tr><tr><td>Learning rate</td><td>0.001</td></tr><tr><td>Epsilon</td><td>0.001</td></tr><tr><td colspan="2">Input features</td></tr><tr><td>pt</td><td>(443), pt (482), pt (561), pt (655), pt (865), VZA, SZA, RAA, θs, RH</td></tr><tr><td colspan="2">Output features</td></tr><tr><td>R<sub>rs</sub></td><td>(443), R<sub>rs</sub> (482), R<sub>rs</sub> (561), R<sub>rs</sub> (655)</td></tr></tbody></table>

<table><caption>**TABLE 3** | PAM-derived aerosol models. G: global, R: regional, Eur: Europe, Amer: America, Aust: Australia, Afr: Africa, N<sub>clus</sub>: number of records per cluster, CS: Cluster Size, CSG: Cluster Structure Goodness.</caption><thead><tr><th>Type</th><th>Sub-type</th><th>Domain</th><th>Continent</th><th>CS(%)</th><th>CSG(%)</th><th>N<sub>clus</sub></th></tr></thead><tbody><tr><td>DU</td><td>1</td><td>R</td><td>Eur, Amer, Afr, Asia, Aust</td><td>3.6</td><td>30</td><td>247</td></tr><tr><td>DU</td><td>2</td><td>R</td><td>Eur, Amer, Afr, Asia, Aust</td><td>3.8</td><td>14</td><td>260</td></tr><tr><td>DU</td><td>3</td><td>R</td><td>Eur, Amer, Afr, Asia, Aust</td><td>3.4</td><td>24</td><td>229</td></tr><tr><td>DU</td><td>4</td><td>R</td><td>Eur, Amer, Afr, Asia</td><td>4.7</td><td>8</td><td>319</td></tr><tr><td>DU</td><td>5</td><td>R</td><td>Eur, Afr, Asia</td><td>2.2</td><td>10</td><td>151</td></tr><tr><td>DU</td><td>6</td><td>R</td><td>Eur, Amer, Afr, Asia, Aust</td><td>3.3</td><td>11</td><td>224</td></tr><tr><td>DU</td><td>7</td><td>R</td><td>Eur, Amer, Afr, Asia, Aust</td><td>1.2</td><td>27</td><td>81</td></tr><tr><td>DU</td><td>8</td><td>R</td><td>Afr, Asia</td><td>4.1</td><td>20</td><td>279</td></tr><tr><td>DU</td><td>9</td><td>R</td><td>Eur, Amer, Afr, Asia, Aust</td><td>5.5</td><td>29</td><td>375</td></tr><tr><td>MAR</td><td>10</td><td>R</td><td>Eur, Amer, Afr, Asia, Aust</td><td>4.4</td><td>19</td><td>301</td></tr><tr><td>MIX</td><td>11</td><td>G</td><td>Eur, Amer, Afr, Asia</td><td>4.5</td><td>18</td><td>306</td></tr><tr><td>MIX</td><td>12</td><td>G</td><td>Eur, Amer, Afr, Asia, Aust</td><td>3.0</td><td>29</td><td>206</td></tr><tr><td>MIX</td><td>13</td><td>R</td><td>Eur, Afr, Asia</td><td>2.0</td><td>9</td><td>139</td></tr><tr><td>MIX</td><td>14</td><td>R</td><td>Eur, Amer, Afr, Asia</td><td>4.6</td><td>24</td><td>310</td></tr><tr><td>MIX</td><td>15</td><td>R</td><td>Eur, Amer, Afr, Asia, Aust</td><td>3.5</td><td>21</td><td>238</td></tr><tr><td>U/I</td><td>16</td><td>G</td><td>Eur, Amer, Afr, Asia</td><td>3.6</td><td>21</td><td>243</td></tr><tr><td>U/I</td><td>17</td><td>R</td><td>Eur, Amer, Afr, Asia, Aust</td><td>3.4</td><td>11</td><td>233</td></tr><tr><td>U/I</td><td>18</td><td>R</td><td>Eur, Amer, Afr, Asia, Aust</td><td>1.5</td><td>29</td><td>101</td></tr><tr><td>U/I</td><td>19</td><td>R</td><td>Eur, Amer, Afr, Asia</td><td>3.6</td><td>30</td><td>242</td></tr><tr><td>U/I</td><td>20</td><td>R</td><td>Eur, Amer, Afr, Asia</td><td>4.1</td><td>27</td><td>279</td></tr><tr><td>U/I</td><td>21</td><td>R</td><td>Eur, Amer, Afr, Asia</td><td>3.8</td><td>12</td><td>259</td></tr><tr><td>BB</td><td>22</td><td>G</td><td>Eur, Amer, Afr, Asia</td><td>4.0</td><td>18</td><td>272</td></tr><tr><td>BB</td><td>23</td><td>G</td><td>Eur, Amer, Afr, Asia</td><td>2.5</td><td>16</td><td>173</td></tr><tr><td>BB</td><td>24</td><td>G</td><td>Eur, Amer, Afr, Asia</td><td>4.1</td><td>16</td><td>279</td></tr><tr><td>BB</td><td>25</td><td>G</td><td>Eur, Amer, Afr, Asia</td><td>4.3</td><td>20</td><td>290</td></tr><tr><td>BB</td><td>26</td><td>R</td><td>Eur, Amer, Afr, Asia</td><td>6.1</td><td>18</td><td>417</td></tr><tr><td>BB</td><td>27</td><td>G</td><td>Eur, Amer, Asia</td><td>4.5</td><td>10</td><td>303</td></tr></tbody></table>

(Omar et al., 2005; Giles et al., 2012; Russell et al., 2014). This aerosol classification was also critical for validating and predicting the monthly and annual composition of dominant aerosols, as described in **section 3.2**. A few representative AERONET sites assigned to different clusters are illustrated in **Figure 4**.

Here, one aerosol type can be grouped into more than one PAM-based cluster, thus, PAM models can be defined for aerosol sub-types that characterize a new mixed type. The PAM models in **Table 3** are categorized according to the mean values of fine mode C<sub>v</sub> fraction [i.e., FC<sub>v</sub><sup>fine</sup> = C<sub>v</sub><sup>fine</sup>/(C<sub>v</sub><sup>fine</sup> + C<sub>v</sub><sup>coarse</sup>)] per aerosol type. Therefore, aerosol sub-types 1 through 9 correspond to DU (C<sub>v</sub> fraction as percentage = 10.2–25%), 10 to MAR (43%), 11–15 to MIX (27.4–59.4%), 16–21 to U/I (34.6–67.5%) and

22–27 to BB (48–65.9%). The CS values among sub-types were highly variable (1.2–6.1%), but the average value per type was comparable between the basic five aerosol sources used here (~3%). The CSG values varied mainly between 8 and 30%. All PAM models had a designated “regional” distribution (i.e., not all continents had measurements); however, the BB aerosol type tended to be represented by a “global” distribution (i.e., AERONET records in all continents but sub-Arctic locations).

In general, most AERONET measurements applied in the PAM classification (95%) were concentrated in the Northern Hemisphere and corresponded to spring-summer months (63.3%) (**Table 4**). Despite this latitudinal asymmetry, BB, MAR, MIX, and U/I aerosol types tended to be more
<table><caption>**TABLE 2** | Properties of MDNs used for developing our experimental AC processors.</caption><thead><tr><th colspan="2">Network hyperparameters</th></tr></thead><tbody><tr><td>Training Iterations</td><td>5,000</td></tr><tr><td>Number of Gaussians</td><td>5</td></tr><tr><td>Neurons per hidden layer</td><td>30</td></tr><tr><td>Hidden layers</td><td>5</td></tr><tr><td>Learning rate</td><td>0.001</td></tr><tr><td>Epsilon</td><td>0.001</td></tr><tr><td colspan="2">Input features</td></tr><tr><td>pt (443), pt (482), pt (561), pt (655), pt (865), VZA, SZA, RAA, θs, RH</td><td></td></tr><tr><td colspan="2">Output features</td></tr><tr><td>R<sub>rs</sub> (443), R<sub>rs</sub> (482), R<sub>rs</sub> (561), R<sub>rs</sub> (655)</td><td></td></tr></tbody></table>

<table><caption>**TABLE 3** | PAM-derived aerosol models. G: global, R: regional, Eur: Europe, Amer: America, Aust: Australia, Afr: Africa, N<sub>clus</sub>: number of records per cluster, CS: Cluster Size, CSG: Cluster Structure Goodness.</caption><thead><tr><th>Type</th><th>Sub-type</th><th>Domain</th><th>Continent</th><th>CS(%)</th><th>CSG(%)</th><th>N<sub>clus</sub></th></tr></thead><tbody><tr><td>DU</td><td>1</td><td>R</td><td>Eur, Amer, Afr, Asia, Aust</td><td>3.6</td><td>30</td><td>247</td></tr><tr><td>DU</td><td>2</td><td>R</td><td>Eur, Amer, Afr, Asia, Aust</td><td>3.8</td><td>14</td><td>260</td></tr><tr><td>DU</td><td>3</td><td>R</td><td>Eur, Amer, Afr, Asia, Aust</td><td>3.4</td><td>24</td><td>229</td></tr><tr><td>DU</td><td>4</td><td>R</td><td>Eur, Amer, Afr, Asia</td><td>4.7</td><td>8</td><td>319</td></tr><tr><td>DU</td><td>5</td><td>R</td><td>Eur, Afr, Asia</td><td>2.2</td><td>10</td><td>151</td></tr><tr><td>DU</td><td>6</td><td>R</td><td>Eur, Amer, Afr, Asia, Aust</td><td>3.3</td><td>11</td><td>224</td></tr><tr><td>DU</td><td>7</td><td>R</td><td>Eur, Amer, Afr, Asia, Aust</td><td>1.2</td><td>27</td><td>81</td></tr><tr><td>DU</td><td>8</td><td>R</td><td>Afr, Asia</td><td>4.1</td><td>20</td><td>279</td></tr><tr><td>DU</td><td>9</td><td>R</td><td>Eur, Amer, Afr, Asia, Aust</td><td>5.5</td><td>29</td><td>375</td></tr><tr><td>MAR</td><td>10</td><td>R</td><td>Eur, Amer, Afr, Asia, Aust</td><td>4.4</td><td>19</td><td>301</td></tr><tr><td>MIX</td><td>11</td><td>G</td><td>Eur, Amer, Afr, Asia</td><td>4.5</td><td>18</td><td>306</td></tr><tr><td>MIX</td><td>12</td><td>G</td><td>Eur, Amer, Afr, Asia, Aust</td><td>3.0</td><td>29</td><td>206</td></tr><tr><td>MIX</td><td>13</td><td>R</td><td>Eur, Afr, Asia</td><td>2.0</td><td>9</td><td>139</td></tr><tr><td>MIX</td><td>14</td><td>R</td><td>Eur, Amer, Afr, Asia</td><td>4.6</td><td>24</td><td>310</td></tr><tr><td>MIX</td><td>15</td><td>R</td><td>Eur, Amer, Afr, Asia, Aust</td><td>3.5</td><td>21</td><td>238</td></tr><tr><td>U/I</td><td>16</td><td>G</td><td>Eur, Amer, Afr, Asia</td><td>3.6</td><td>21</td><td>243</td></tr><tr><td>U/I</td><td>17</td><td>R</td><td>Eur, Amer, Afr, Asia, Aust</td><td>3.4</td><td>11</td><td>233</td></tr><tr><td>U/I</td><td>18</td><td>R</td><td>Eur, Amer, Afr, Asia, Aust</td><td>1.5</td><td>29</td><td>101</td></tr><tr><td>U/I</td><td>19</td><td>R</td><td>Eur, Amer, Afr, Asia</td><td>3.6</td><td>30</td><td>242</td></tr><tr><td>U/I</td><td>20</td><td>R</td><td>Eur, Amer, Afr, Asia</td><td>4.1</td><td>27</td><td>279</td></tr><tr><td>U/I</td><td>21</td><td>R</td><td>Eur, Amer, Afr, Asia</td><td>3.8</td><td>12</td><td>259</td></tr><tr><td>BB</td><td>22</td><td>G</td><td>Eur, Amer, Afr, Asia</td><td>4.0</td><td>18</td><td>272</td></tr><tr><td>BB</td><td>23</td><td>G</td><td>Eur, Amer, Afr, Asia</td><td>2.5</td><td>16</td><td>173</td></tr><tr><td>BB</td><td>24</td><td>G</td><td>Eur, Amer, Afr, Asia</td><td>4.1</td><td>16</td><td>279</td></tr><tr><td>BB</td><td>25</td><td>G</td><td>Eur, Amer, Afr, Asia</td><td>4.3</td><td>20</td><td>290</td></tr><tr><td>BB</td><td>26</td><td>R</td><td>Eur, Amer, Afr, Asia</td><td>6.1</td><td>18</td><td>417</td></tr><tr><td>BB</td><td>27</td><td>G</td><td>Eur, Amer, Asia</td><td>4.5</td><td>10</td><td>303</td></tr></tbody></table>

(Omar et al., 2005; Giles et al., 2012; Russell et al., 2014). This aerosol classification was also critical for validating and predicting the monthly and annual composition of dominant aerosols, as described in **section 3.2**. A few representative AERONET sites assigned to different clusters are illustrated in **Figure 4**.

Here, one aerosol type can be grouped into more than one PAM-based cluster, thus, PAM models can be defined for aerosol sub-types that characterize a new mixed type. The PAM models in **Table 3** are categorized according to the mean values of fine mode $C_v$ fraction [i.e., $FC_v^{fine} = C_v^{fine}/(C_v^{fine} + C_v^{coarse})$] per aerosol type. Therefore, aerosol sub-types 1 through 9 correspond to DU ($C_v$ fraction as percentage = 10.2–25%), 10 to MAR (43%), 11–15 to MIX (27.4–59.4%), 16–21 to U/I (34.6–67.5%) and

22–27 to BB (48–65.9%). The CS values among sub-types were highly variable (1.2–6.1%), but the average value per type was comparable between the basic five aerosol sources used here (~3%). The CSG values varied mainly between 8 and 30%. All PAM models had a designated “regional” distribution (i.e., not all continents had measurements); however, the BB aerosol type tended to be represented by a “global” distribution (i.e., AERONET records in all continents but sub-Arctic locations).

In general, most AERONET measurements applied in the PAM classification (95%) were concentrated in the Northern Hemisphere and corresponded to spring-summer months (63.3%) (**Table 4**). Despite this latitudinal asymmetry, BB, MAR, MIX, and U/I aerosol types tended to be more
<table><caption>**TABLE 4** | Seasonal aerosol-type distributions. Values represent the total percentage (i.e., Northern and Southern Hemispheres) for each aerosol category. Each column sums up to 100%.</caption><thead><tr><th></th><th>Winter</th><th>Spring</th><th>Summer</th><th>Fall</th></tr></thead><tbody><tr><th>DU</th><td>35.1</td><td>43.1</td><td>26.6</td><td>28.3</td></tr><tr><th>MAR</th><td>5.5</td><td>5.4</td><td>3.4</td><td>5.3</td></tr><tr><th>MIX</th><td>22.3</td><td>18.2</td><td>16.0</td><td>17.7</td></tr><tr><th>U/I</th><td>17.6</td><td>16.4</td><td>22.8</td><td>20.4</td></tr><tr><th>BB</th><td>19.5</td><td>16.9</td><td>31.2</td><td>28.3</td></tr></tbody></table>

prevalent over sites south of the equator. This spatial pattern was reversed for DU, possibly due to both the lack of significant dust sources (e.g., compared to Saharan and Gobi Deserts) and the location of AERONET sites near dust regions in the Southern Hemisphere.

In general, DU models were dominant in all months, representing 32% of aerosol records with the highest contribution during spring (43.1%). The AERONET measurements with biomass burning classification were the second most dominant aerosol type during an annual cycle with the highest contributions during summer (31.2%), followed closely by the fall period (28.3%) (Table 4). Marine-dominated aerosols only account for less than 6% in all seasons, possibly due to the limited availability of inverted AERONET aerosol properties for AOD(0.44) < 0.4 across most of these locations. The highest seasonal contribution associated with the MIX aerosol type was 22.3% and occurred during winter. Also, the influence of MIX aerosols decreased during warmer months and was characterized by a weak inter-seasonal variability. Like BB, the U/I aerosol type had the highest seasonal contribution of (22.8%) during summer followed by the fall period (20.4%).

**Table 5** summarizes the mean values of SSA(0.44), AE<sub>e</sub> and AE<sub>a</sub> for each aerosol class derived from PAM and the mean values of other relevant aerosol properties (e.g., g, r<sub>n</sub>, and n<sub>r</sub>) correlated with our input classification parameters. These values are within the expected range for aerosol types described in the literature (Dubovik et al., 2002a; Giles et al., 2012; Russell et al., 2014; Roger et al., 2022). The standard deviation of each aerosol property is not shown as the intra-cluster variability can be inferred from CSG values computed in Table 3. In general, the magnitude of AE<sub>e</sub> decreased (increased) at lower (higher) aerosol sub-type numbers (e.g., 0.27 to 0.70 for DU sub-types 1–9) (e.g., 1.46 to 1.88 for BB sub-types 22–27) due to smaller (larger) contributions of fine aerosol particles. Mean values of AE<sub>a</sub> for DU ranged from 1.5 to 2.5, with values closer to 1 typically associated with the remaining aerosol types. On average, fine and coarse mode particles had the smallest radius for DU sub-types (0.04 and 0.69 µm, respectively). In general, σ<sub>n</sub><sup>fine</sup> (σ<sub>n</sub><sup>coarse</sup>) decreased (increased) in aerosol types with higher FCV<sup>fine</sup> values (e.g., σ<sub>n</sub><sup>fine</sup> up to 0.63 and 0.45 in DU and BB, respectively). The real (imaginary) part of the refractive index for DU-like models was commonly high (low) with respect to the other aerosol types (e.g., maximum n<sub>r</sub>(0.44) and minimum n<sub>i</sub>(0.44) for DU sub-types were 1.498 and 0.001, respectively). The variability for n<sub>r</sub> (n<sub>i</sub>) in

<figure><img src="image_5.png" /><figcaption>**FIGURE 4** | Examples of PAM-derived aerosol types in AERONET sites. Sub-types 1–10. Examples of PAM-derived aerosol types in AERONET sites. Sub-types 11–21. Examples of PAM-derived aerosol types in AERONET sites. Sub-types 22–27.</figcaption></figure>
<table><caption>**TABLE 4** | Seasonal aerosol-type distributions. Values represent the total percentage (i.e., Northern and Southern Hemispheres) for each aerosol category. Each column sums up to 100%.</caption><thead><tr><th></th><th>Winter</th><th>Spring</th><th>Summer</th><th>Fall</th></tr></thead><tbody><tr><th>DU</th><td>35.1</td><td>43.1</td><td>26.6</td><td>28.3</td></tr><tr><th>MAR</th><td>5.5</td><td>5.4</td><td>3.4</td><td>5.3</td></tr><tr><th>MIX</th><td>22.3</td><td>18.2</td><td>16.0</td><td>17.7</td></tr><tr><th>U/I</th><td>17.6</td><td>16.4</td><td>22.8</td><td>20.4</td></tr><tr><th>BB</th><td>19.5</td><td>16.9</td><td>31.2</td><td>28.3</td></tr></tbody></table>

prevalent over sites south of the equator. This spatial pattern was reversed for DU, possibly due to both the lack of significant dust sources (e.g., compared to Saharan and Gobi Deserts) and the location of AERONET sites near dust regions in the Southern Hemisphere.

In general, DU models were dominant in all months, representing 32% of aerosol records with the highest contribution during spring (43.1%). The AERONET measurements with biomass burning classification were the second most dominant aerosol type during an annual cycle with the highest contributions during summer (31.2%), followed closely by the fall period (28.3%) (Table 4). Marine-dominated aerosols only account for less than 6% in all seasons, possibly due to the limited availability of inverted AERONET aerosol properties for AOD(0.44) < 0.4 across most of these locations. The highest seasonal contribution associated with the MIX aerosol type was 22.3% and occurred during winter. Also, the influence of MIX aerosols decreased during warmer months and was characterized by a weak inter-seasonal variability. Like BB, the U/I aerosol type had the highest seasonal contribution of (22.8%) during summer followed by the fall period (20.4%).

**Table 5** summarizes the mean values of SSA(0.44), AE<sub>e</sub> and AE<sub>a</sub> for each aerosol class derived from PAM and the mean values of other relevant aerosol properties (e.g., g, r<sub>n</sub>, and n<sub>r</sub>) correlated with our input classification parameters. These values are within the expected range for aerosol types described in the literature (Dubovik et al., 2002a; Giles et al., 2012; Russell et al., 2014; Roger et al., 2022). The standard deviation of each aerosol property is not shown as the intra-cluster variability can be inferred from CSG values computed in Table 3. In general, the magnitude of AE<sub>e</sub> decreased (increased) at lower (higher) aerosol sub-type numbers (e.g., 0.27 to 0.70 for DU sub-types 1–9) (e.g., 1.46 to 1.88 for BB sub-types 22–27) due to smaller (larger) contributions of fine aerosol particles. Mean values of AE<sub>a</sub> for DU ranged from 1.5 to 2.5, with values closer to 1 typically associated with the remaining aerosol types. On average, fine and coarse mode particles had the smallest radius for DU sub-types (0.04 and 0.69 µm, respectively). In general, σ<sub>n</sub><sup>fine</sup> (σ<sub>n</sub><sup>coarse</sup>) decreased (increased) in aerosol types with higher FCV<sup>fine</sup> values (e.g., σ<sub>n</sub><sup>fine</sup> up to 0.63 and 0.45 in DU and BB, respectively). The real (imaginary) part of the refractive index for DU-like models was commonly high (low) with respect to the other aerosol types (e.g., maximum n<sub>r</sub>(0.44) and minimum n<sub>i</sub>(0.44) for DU sub-types were 1.498 and 0.001, respectively). The variability for n<sub>r</sub> (n<sub>i</sub>) in

<figure><img src="image_6.png" /><figcaption>**FIGURE 4** | Examples of PAM-derived aerosol types in AERONET sites. Sub-types 1–10. Examples of PAM-derived aerosol types in AERONET sites. Sub-types 11–21. Examples of PAM-derived aerosol types in AERONET sites. Sub-types 22–27.</figcaption></figure>
<table><caption>**TABLE 5** | Summary of aerosol properties per cluster. Each value corresponds to the arithmetic mean. Units of properties are described in **Table 1**. Sub-types 1–9 are DU, 10 is MAR, 11–15 are MIX, 16–21 are U/I, and 22–27 are BB.</caption><thead><tr><th>Sub-type</th><th>SSA(0.44)</th><th>AE<sub>a</sub></th><th>AE<sub>a</sub></th><th>FC<sub>v</sub><sup>fine</sup></th><th>r<sub>n</sub><sup>fine</sup></th><th>σ<sub>n</sub><sup>fine</sup></th><th>r<sub>n</sub><sup>coarse</sup></th><th>σ<sub>n</sub><sup>coarse</sup></th><th>n<sub>r</sub>(0.44)</th><th>n<sub>f</sub>(0.44)</th><th>n<sub>r</sub>(0.87)</th><th>n<sub>f</sub>(0.87)</th><th>SSA(0.87)</th><th>g(0.44)</th><th>g(0.87)</th></tr></thead><tbody><tr><td>1</td><td>0.93</td><td>0.273</td><td>2.572</td><td>0.11</td><td>0.04</td><td>0.63</td><td>0.73</td><td>0.58</td><td>1.465</td><td>0.003</td><td>1.470</td><td>0.001</td><td>0.98</td><td>0.773</td><td>0.742</td></tr><tr><td>2</td><td>0.86</td><td>0.600</td><td>1.640</td><td>0.17</td><td>0.06</td><td>0.52</td><td>0.72</td><td>0.64</td><td>1.488</td><td>0.009</td><td>1.525</td><td>0.005</td><td>0.93</td><td>0.744</td><td>0.715</td></tr><tr><td>3</td><td>0.93</td><td>0.593</td><td>2.186</td><td>0.17</td><td>0.06</td><td>0.52</td><td>0.70</td><td>0.62</td><td>1.465</td><td>0.004</td><td>1.463</td><td>0.002</td><td>0.97</td><td>0.737</td><td>0.712</td></tr><tr><td>4</td><td>0.88</td><td>0.343</td><td>1.990</td><td>0.11</td><td>0.05</td><td>0.60</td><td>0.75</td><td>0.6</td><td>1.481</td><td>0.006</td><td>1.502</td><td>0.003</td><td>0.96</td><td>0.774</td><td>0.735</td></tr><tr><td>5</td><td>0.82</td><td>0.707</td><td>1.421</td><td>0.19</td><td>0.06</td><td>0.51</td><td>0.70</td><td>0.65</td><td>1.498</td><td>0.016</td><td>1.532</td><td>0.008</td><td>0.89</td><td>0.743</td><td>0.711</td></tr><tr><td>6</td><td>0.90</td><td>0.963</td><td>1.452</td><td>0.25</td><td>0.08</td><td>0.47</td><td>0.69</td><td>0.66</td><td>1.496</td><td>0.008</td><td>1.510</td><td>0.006</td><td>0.92</td><td>0.718</td><td>0.681</td></tr><tr><td>7</td><td>0.97</td><td>0.370</td><td>1.637</td><td>0.12</td><td>0.06</td><td>0.58</td><td>0.73</td><td>0.58</td><td>1.468</td><td>0.001</td><td>1.478</td><td>0.001</td><td>0.98</td><td>0.754</td><td>0.731</td></tr><tr><td>8</td><td>0.91</td><td>0.687</td><td>1.772</td><td>0.19</td><td>0.07</td><td>0.49</td><td>0.70</td><td>0.63</td><td>1.471</td><td>0.005</td><td>1.503</td><td>0.003</td><td>0.95</td><td>0.732</td><td>0.706</td></tr><tr><td>9</td><td>0.90</td><td>0.270</td><td>2.463</td><td>0.10</td><td>0.04</td><td>0.63</td><td>0.75</td><td>0.59</td><td>1.473</td><td>0.004</td><td>1.482</td><td>0.002</td><td>0.97</td><td>0.779</td><td>0.741</td></tr><tr><td>10</td><td>0.93</td><td>1.289</td><td>1.141</td><td>0.43</td><td>0.10</td><td>0.46</td><td>0.76</td><td>0.65</td><td>1.466</td><td>0.007</td><td>1.468</td><td>0.007</td><td>0.91</td><td>0.712</td><td>0.638</td></tr><tr><td>11</td><td>0.96</td><td>1.384</td><td>0.975</td><td>0.51</td><td>0.10</td><td>0.46</td><td>0.77</td><td>0.65</td><td>1.450</td><td>0.005</td><td>1.447</td><td>0.006</td><td>0.93</td><td>0.715</td><td>0.628</td></tr><tr><td>12</td><td>0.85</td><td>1.326</td><td>1.337</td><td>0.42</td><td>0.08</td><td>0.47</td><td>0.78</td><td>0.67</td><td>1.493</td><td>0.023</td><td>1.514</td><td>0.017</td><td>0.82</td><td>0.693</td><td>0.620</td></tr><tr><td>13</td><td>0.97</td><td>0.859</td><td>1.237</td><td>0.27</td><td>0.08</td><td>0.49</td><td>0.71</td><td>0.63</td><td>1.447</td><td>0.002</td><td>1.459</td><td>0.002</td><td>0.97</td><td>0.730</td><td>0.694</td></tr><tr><td>14</td><td>0.97</td><td>1.615</td><td>0.980</td><td>0.59</td><td>0.10</td><td>0.43</td><td>0.78</td><td>0.66</td><td>1.441</td><td>0.004</td><td>1.439</td><td>0.005</td><td>0.94</td><td>0.709</td><td>0.584</td></tr><tr><td>15</td><td>0.94</td><td>1.031</td><td>1.389</td><td>0.32</td><td>0.09</td><td>0.47</td><td>0.70</td><td>0.65</td><td>1.462</td><td>0.005</td><td>1.468</td><td>0.004</td><td>0.94</td><td>0.721</td><td>0.674</td></tr><tr><td>16</td><td>0.88</td><td>1.184</td><td>1.250</td><td>0.35</td><td>0.08</td><td>0.47</td><td>0.75</td><td>0.67</td><td>1.502</td><td>0.014</td><td>1.517</td><td>0.012</td><td>0.87</td><td>0.703</td><td>0.649</td></tr><tr><td>17</td><td>0.98</td><td>1.238</td><td>0.967</td><td>0.44</td><td>0.10</td><td>0.46</td><td>0.76</td><td>0.63</td><td>1.446</td><td>0.002</td><td>1.446</td><td>0.002</td><td>0.97</td><td>0.726</td><td>0.653</td></tr><tr><td>18</td><td>0.80</td><td>1.395</td><td>1.339</td><td>0.45</td><td>0.08</td><td>0.47</td><td>0.81</td><td>0.67</td><td>1.496</td><td>0.033</td><td>1.520</td><td>0.026</td><td>0.76</td><td>0.680</td><td>0.604</td></tr><tr><td>19</td><td>0.99</td><td>1.821</td><td>1.126</td><td>0.68</td><td>0.11</td><td>0.41</td><td>0.86</td><td>0.64</td><td>1.417</td><td>0.002</td><td>1.425</td><td>0.002</td><td>0.97</td><td>0.705</td><td>0.544</td></tr><tr><td>20</td><td>0.99</td><td>1.289</td><td>1.062</td><td>0.58</td><td>0.10</td><td>0.45</td><td>0.83</td><td>0.63</td><td>1.421</td><td>0.002</td><td>1.429</td><td>0.002</td><td>0.98</td><td>0.719</td><td>0.603</td></tr><tr><td>21</td><td>0.95</td><td>1.924</td><td>1.077</td><td>0.65</td><td>0.10</td><td>0.40</td><td>0.85</td><td>0.66</td><td>1.471</td><td>0.007</td><td>1.453</td><td>0.008</td><td>0.90</td><td>0.686</td><td>0.525</td></tr><tr><td>22</td><td>0.93</td><td>1.819</td><td>1.094</td><td>0.61</td><td>0.10</td><td>0.41</td><td>0.84</td><td>0.66</td><td>1.481</td><td>0.011</td><td>1.463</td><td>0.012</td><td>0.86</td><td>0.683</td><td>0.540</td></tr><tr><td>23</td><td>0.86</td><td>1.731</td><td>1.294</td><td>0.58</td><td>0.09</td><td>0.42</td><td>0.83</td><td>0.68</td><td>1.502</td><td>0.025</td><td>1.513</td><td>0.022</td><td>0.78</td><td>0.659</td><td>0.532</td></tr><tr><td>24</td><td>0.90</td><td>1.651</td><td>1.192</td><td>0.55</td><td>0.09</td><td>0.43</td><td>0.82</td><td>0.67</td><td>1.473</td><td>0.015</td><td>1.470</td><td>0.015</td><td>0.83</td><td>0.681</td><td>0.569</td></tr><tr><td>25</td><td>0.92</td><td>1.456</td><td>1.119</td><td>0.48</td><td>0.09</td><td>0.45</td><td>0.76</td><td>0.66</td><td>1.471</td><td>0.011</td><td>1.467</td><td>0.011</td><td>0.88</td><td>0.703</td><td>0.614</td></tr><tr><td>26</td><td>0.94</td><td>1.631</td><td>0.996</td><td>0.57</td><td>0.10</td><td>0.43</td><td>0.82</td><td>0.65</td><td>1.457</td><td>0.008</td><td>1.445</td><td>0.008</td><td>0.90</td><td>0.704</td><td>0.583</td></tr><tr><td>27</td><td>0.97</td><td>1.879</td><td>1.114</td><td>0.66</td><td>0.11</td><td>0.40</td><td>0.87</td><td>0.65</td><td>1.445</td><td>0.004</td><td>1.435</td><td>0.004</td><td>0.94</td><td>0.699</td><td>0.538</td></tr></tbody></table>

the blue or NIR spectral range was relatively small (large) with values < 5% and >90%, respectively. As expected, the asymmetry parameter tended to decrease as FC<sub>v</sub><sup>fine</sup> increases (e.g., the mean g(0.44) for DU and BB types was 0.750 and 0.688, respectively).

The differences between PAM-derived aerosol types are also present when angular and spectral changes of aerosol properties are considered. To emphasize these variations, selected models represented in the four aerosol types are depicted in **Figure 5**. In general, DU-MIX (BB-U/I) sub-types had the lowest (highest) backscattering contributions as inferred from their normalized phase (P<sub>a</sub>) functions (i.e., θ<sub>s</sub> = 90° to 160°) (**Figure 5B**) with the MAR models having intermediate values. Note that several curves may be depicted with the same color when the main aerosol type is composed of more than one PAM model (i.e., sub-type).

In the forward-scattering direction (i.e., θ<sub>s</sub> = 5 to 25°) (**Figure 5A**), there was not a clear differentiation between DU-MIX and BB-U/I sub-types. Likewise, aerosols mainly composed of sea salt were characterized by a normalized phase function falling in the middle of the envelope composed of all aerosol models. The spectral behavior of SSA was a good indicator to discriminate aerosols with different compositions as three general patterns are observed: P1. aerosols with SSA increasing with λ (e.g., DU), P2. aerosols with SSA decreasing with λ (e.g., BB and MAR), and P3. aerosols presenting P1 and P2 (e.g., U/I and MIX) (**Figure 5C**). Spectral curves of g are depicted in **Figure 5D** and suggest a direct relationship between the spectral slope and FC<sub>v</sub><sup>fine</sup>. More specifically, the spectral slope of g between 0.44 and 1.02 µm were high (low) for fine-mode

(coarse-mode) dominated aerosol assemblages (i.e., BB-U/I and DU sub-types, respectively). Likewise, aerosols with a substantial contribution of fine and coarse particles (i.e., MIX and MAR) presented intermediate values of g spectral slopes. The complementary components of the PAM models are readily inferred by comparing these spectra with the A2010 models shown in **Figure 1**.

## Model Comparisons With Literature

The calculation of dominant clusters over AERONET sites reported in the literature for different aerosol types was used to evaluate the PAM performance for identifying aerosol reference sources. This verification was based on annually dominant aerosol categories. Thus, AERONET records corresponding to different periods of the year were merged for each annual estimate. This comparison technique had varying degrees of success as the proportion of correct matches for DU, BB, U/I, MIX, and MAR were highly variable (i.e., 67, 33, 32, 15, and 10% success rate, respectively). In general, AERONET sites with correct matches for DU (e.g., Banizoumbou, Cape Verde, Ouagadougou) were previously proposed by more than one study (see **Supplementary Appendix Table A3**). The sites with accurate classification of BB aerosols (e.g., Mongu, Bonanza Creek, Cuiaba) correspond with aerosol sub-type clusters 23 to 25. The PAM-based classification of the GSFC site suggested aerosols mainly originated from urban/industrial activities, which agreed with reference aerosol sources cited (**Supplementary Appendix Table A3**). Most sites with correct matches in U/I
<table><caption>TABLE 5 | Summary of aerosol properties per cluster. Each value corresponds to the arithmetic mean. Units of properties are described in Table 1. Sub-types 1–9 are DU, 10 is MAR, 11–15 are MIX, 16–21 are U/I, and 22–27 are BB.</caption><thead><tr><th>Sub-type</th><th>SSA(0.44)</th><th>AE<sub>a</sub></th><th>AE<sub>a</sub></th><th>FC<sub>v</sub><sup>fine</sup></th><th>r<sub>n</sub><sup>fine</sup></th><th>σ<sub>n</sub><sup>fine</sup></th><th>r<sub>n</sub><sup>coarse</sup></th><th>σ<sub>n</sub><sup>coarse</sup></th><th>n<sub>r</sub>(0.44)</th><th>n<sub>f</sub>(0.44)</th><th>n<sub>r</sub>(0.87)</th><th>n<sub>f</sub>(0.87)</th><th>SSA(0.87)</th><th>g(0.44)</th><th>g(0.87)</th></tr></thead><tbody><tr><td>1</td><td>0.93</td><td>0.273</td><td>2.572</td><td>0.11</td><td>0.04</td><td>0.63</td><td>0.73</td><td>0.58</td><td>1.465</td><td>0.003</td><td>1.470</td><td>0.001</td><td>0.98</td><td>0.773</td><td>0.742</td></tr><tr><td>2</td><td>0.86</td><td>0.600</td><td>1.640</td><td>0.17</td><td>0.06</td><td>0.52</td><td>0.72</td><td>0.64</td><td>1.488</td><td>0.009</td><td>1.525</td><td>0.005</td><td>0.93</td><td>0.744</td><td>0.715</td></tr><tr><td>3</td><td>0.93</td><td>0.593</td><td>2.186</td><td>0.17</td><td>0.06</td><td>0.52</td><td>0.70</td><td>0.62</td><td>1.465</td><td>0.004</td><td>1.463</td><td>0.002</td><td>0.97</td><td>0.737</td><td>0.712</td></tr><tr><td>4</td><td>0.88</td><td>0.343</td><td>1.990</td><td>0.11</td><td>0.05</td><td>0.60</td><td>0.75</td><td>0.6</td><td>1.481</td><td>0.006</td><td>1.502</td><td>0.003</td><td>0.96</td><td>0.774</td><td>0.735</td></tr><tr><td>5</td><td>0.82</td><td>0.707</td><td>1.421</td><td>0.19</td><td>0.06</td><td>0.51</td><td>0.70</td><td>0.65</td><td>1.498</td><td>0.016</td><td>1.532</td><td>0.008</td><td>0.89</td><td>0.743</td><td>0.711</td></tr><tr><td>6</td><td>0.90</td><td>0.963</td><td>1.452</td><td>0.25</td><td>0.08</td><td>0.47</td><td>0.69</td><td>0.66</td><td>1.496</td><td>0.008</td><td>1.510</td><td>0.006</td><td>0.92</td><td>0.718</td><td>0.681</td></tr><tr><td>7</td><td>0.97</td><td>0.370</td><td>1.637</td><td>0.12</td><td>0.06</td><td>0.58</td><td>0.73</td><td>0.58</td><td>1.468</td><td>0.001</td><td>1.478</td><td>0.001</td><td>0.98</td><td>0.754</td><td>0.731</td></tr><tr><td>8</td><td>0.91</td><td>0.687</td><td>1.772</td><td>0.19</td><td>0.07</td><td>0.49</td><td>0.70</td><td>0.63</td><td>1.471</td><td>0.005</td><td>1.503</td><td>0.003</td><td>0.95</td><td>0.732</td><td>0.706</td></tr><tr><td>9</td><td>0.90</td><td>0.270</td><td>2.463</td><td>0.10</td><td>0.04</td><td>0.63</td><td>0.75</td><td>0.59</td><td>1.473</td><td>0.004</td><td>1.482</td><td>0.002</td><td>0.97</td><td>0.779</td><td>0.741</td></tr><tr><td>10</td><td>0.93</td><td>1.289</td><td>1.141</td><td>0.43</td><td>0.10</td><td>0.46</td><td>0.76</td><td>0.65</td><td>1.466</td><td>0.007</td><td>1.468</td><td>0.007</td><td>0.91</td><td>0.712</td><td>0.638</td></tr><tr><td>11</td><td>0.96</td><td>1.384</td><td>0.975</td><td>0.51</td><td>0.10</td><td>0.46</td><td>0.77</td><td>0.65</td><td>1.450</td><td>0.005</td><td>1.447</td><td>0.006</td><td>0.93</td><td>0.715</td><td>0.628</td></tr><tr><td>12</td><td>0.85</td><td>1.326</td><td>1.337</td><td>0.42</td><td>0.08</td><td>0.47</td><td>0.78</td><td>0.67</td><td>1.493</td><td>0.023</td><td>1.514</td><td>0.017</td><td>0.82</td><td>0.693</td><td>0.620</td></tr><tr><td>13</td><td>0.97</td><td>0.859</td><td>1.237</td><td>0.27</td><td>0.08</td><td>0.49</td><td>0.71</td><td>0.63</td><td>1.447</td><td>0.002</td><td>1.459</td><td>0.002</td><td>0.97</td><td>0.730</td><td>0.694</td></tr><tr><td>14</td><td>0.97</td><td>1.615</td><td>0.980</td><td>0.59</td><td>0.10</td><td>0.43</td><td>0.78</td><td>0.66</td><td>1.441</td><td>0.004</td><td>1.439</td><td>0.005</td><td>0.94</td><td>0.709</td><td>0.584</td></tr><tr><td>15</td><td>0.94</td><td>1.031</td><td>1.389</td><td>0.32</td><td>0.09</td><td>0.47</td><td>0.70</td><td>0.65</td><td>1.462</td><td>0.005</td><td>1.468</td><td>0.004</td><td>0.94</td><td>0.721</td><td>0.674</td></tr><tr><td>16</td><td>0.88</td><td>1.184</td><td>1.250</td><td>0.35</td><td>0.08</td><td>0.47</td><td>0.75</td><td>0.67</td><td>1.502</td><td>0.014</td><td>1.517</td><td>0.012</td><td>0.87</td><td>0.703</td><td>0.649</td></tr><tr><td>17</td><td>0.98</td><td>1.238</td><td>0.967</td><td>0.44</td><td>0.10</td><td>0.46</td><td>0.76</td><td>0.63</td><td>1.446</td><td>0.002</td><td>1.446</td><td>0.002</td><td>0.97</td><td>0.726</td><td>0.653</td></tr><tr><td>18</td><td>0.80</td><td>1.395</td><td>1.339</td><td>0.45</td><td>0.08</td><td>0.47</td><td>0.81</td><td>0.67</td><td>1.496</td><td>0.033</td><td>1.520</td><td>0.026</td><td>0.76</td><td>0.680</td><td>0.604</td></tr><tr><td>19</td><td>0.99</td><td>1.821</td><td>1.126</td><td>0.68</td><td>0.11</td><td>0.41</td><td>0.86</td><td>0.64</td><td>1.417</td><td>0.002</td><td>1.425</td><td>0.002</td><td>0.97</td><td>0.705</td><td>0.544</td></tr><tr><td>20</td><td>0.99</td><td>1.289</td><td>1.062</td><td>0.58</td><td>0.10</td><td>0.45</td><td>0.83</td><td>0.63</td><td>1.421</td><td>0.002</td><td>1.429</td><td>0.002</td><td>0.98</td><td>0.719</td><td>0.603</td></tr><tr><td>21</td><td>0.95</td><td>1.924</td><td>1.077</td><td>0.65</td><td>0.10</td><td>0.40</td><td>0.85</td><td>0.66</td><td>1.471</td><td>0.007</td><td>1.453</td><td>0.008</td><td>0.90</td><td>0.686</td><td>0.525</td></tr><tr><td>22</td><td>0.93</td><td>1.819</td><td>1.094</td><td>0.61</td><td>0.10</td><td>0.41</td><td>0.84</td><td>0.66</td><td>1.481</td><td>0.011</td><td>1.463</td><td>0.012</td><td>0.86</td><td>0.683</td><td>0.540</td></tr><tr><td>23</td><td>0.86</td><td>1.731</td><td>1.294</td><td>0.58</td><td>0.09</td><td>0.42</td><td>0.83</td><td>0.68</td><td>1.502</td><td>0.025</td><td>1.513</td><td>0.022</td><td>0.78</td><td>0.659</td><td>0.532</td></tr><tr><td>24</td><td>0.90</td><td>1.651</td><td>1.192</td><td>0.55</td><td>0.09</td><td>0.43</td><td>0.82</td><td>0.67</td><td>1.473</td><td>0.015</td><td>1.470</td><td>0.015</td><td>0.83</td><td>0.681</td><td>0.569</td></tr><tr><td>25</td><td>0.92</td><td>1.456</td><td>1.119</td><td>0.48</td><td>0.09</td><td>0.45</td><td>0.76</td><td>0.66</td><td>1.471</td><td>0.011</td><td>1.467</td><td>0.011</td><td>0.88</td><td>0.703</td><td>0.614</td></tr><tr><td>26</td><td>0.94</td><td>1.631</td><td>0.996</td><td>0.57</td><td>0.10</td><td>0.43</td><td>0.82</td><td>0.65</td><td>1.457</td><td>0.008</td><td>1.445</td><td>0.008</td><td>0.90</td><td>0.704</td><td>0.583</td></tr><tr><td>27</td><td>0.97</td><td>1.879</td><td>1.114</td><td>0.66</td><td>0.11</td><td>0.40</td><td>0.87</td><td>0.65</td><td>1.445</td><td>0.004</td><td>1.435</td><td>0.004</td><td>0.94</td><td>0.699</td><td>0.538</td></tr></tbody></table>

the blue or NIR spectral range was relatively small (large) with values < 5% and >90%, respectively. As expected, the asymmetry parameter tended to decrease as FC<sub>v</sub><sup>fine</sup> increases (e.g., the mean g(0.44) for DU and BB types was 0.750 and 0.688, respectively).

The differences between PAM-derived aerosol types are also present when angular and spectral changes of aerosol properties are considered. To emphasize these variations, selected models represented in the four aerosol types are depicted in Figure 5. In general, DU-MIX (BB-U/I) sub-types had the lowest (highest) backscattering contributions as inferred from their normalized phase (P<sub>a</sub>) functions (i.e., θ<sub>s</sub> = 90° to 160°) (Figure 5B) with the MAR models having intermediate values. Note that several curves may be depicted with the same color when the main aerosol type is composed of more than one PAM model (i.e., sub-type).

In the forward-scattering direction (i.e., θ<sub>s</sub> = 5 to 25°) (Figure 5A), there was not a clear differentiation between DU-MIX and BB-U/I sub-types. Likewise, aerosols mainly composed of sea salt were characterized by a normalized phase function falling in the middle of the envelope composed of all aerosol models. The spectral behavior of SSA was a good indicator to discriminate aerosols with different compositions as three general patterns are observed: P1. aerosols with SSA increasing with λ (e.g., DU), P2. aerosols with SSA decreasing with λ (e.g., BB and MAR), and P3. aerosols presenting P1 and P2 (e.g., U/I and MIX) (Figure 5C). Spectral curves of g are depicted in Figure 5D and suggest a direct relationship between the spectral slope and FC<sub>v</sub><sup>fine</sup>. More specifically, the spectral slope of g between 0.44 and 1.02 µm were high (low) for fine-mode

(coarse-mode) dominated aerosol assemblages (i.e., BB-U/I and DU sub-types, respectively). Likewise, aerosols with a substantial contribution of fine and coarse particles (i.e., MIX and MAR) presented intermediate values of g spectral slopes. The complementary components of the PAM models are readily inferred by comparing these spectra with the A2010 models shown in Figure 1.

## Model Comparisons With Literature

The calculation of dominant clusters over AERONET sites reported in the literature for different aerosol types was used to evaluate the PAM performance for identifying aerosol reference sources. This verification was based on annually dominant aerosol categories. Thus, AERONET records corresponding to different periods of the year were merged for each annual estimate. This comparison technique had varying degrees of success as the proportion of correct matches for DU, BB, U/I, MIX, and MAR were highly variable (i.e., 67, 33, 32, 15, and 10% success rate, respectively). In general, AERONET sites with correct matches for DU (e.g., Banizoumbou, Cape Verde, Ouagadougou) were previously proposed by more than one study (see Supplementary Appendix Table A3). The sites with accurate classification of BB aerosols (e.g., Mongu, Bonanza Creek, Cuiaba) correspond with aerosol sub-type clusters 23 to 25. The PAM-based classification of the GSFC site suggested aerosols mainly originated from urban/industrial activities, which agreed with reference aerosol sources cited (Supplementary Appendix Table A3). Most sites with correct matches in U/I
<figure><img src="image_6.png" /><figcaption>**FIGURE 5** | Spectral variation of aerosol properties for selected PAM-derived clusters and associated types. Arithmetic average of normalized phase function at scattering angle $\theta_s = 0^{\circ}$) for forward panel (A) and backscattering angles [panel (B)], SSA [panel (C)], and g [panel (D)]. Curves with the same color symbol represent different aerosol models of the same type (Table 5). The numbers next to aerosol types represent PAM-derived subtypes (Table 5).</figcaption></figure>

coincided with locations studied by Omar et al. (2005) in **Supplementary Appendix Table A3** (e.g., Venise, UCLA, Mexico City) and corresponded to sub-types 20 and 21. The identification of MIX-dominated sites based on PAM was not satisfactory as PAM predictions were only correct in two reference AERONET locations (Kanpur and Howland).

The AERONET sites labeled as DU by the PAM classifier coincided geographically with deserts (**Supplementary Appendix Table A4** and **Figure 4**) (e.g., Mojave, Sonoran). To date, no AERONET sites exist in the Gobi Desert. However, its influence was manifested in western China. Predominance of dust along with observations made in Florida and the Caribbean Sea during summer was likely related to trans-Atlantic transport of dust originated in Sahara and Sub-Sahara (Sahel) deserts of Western Africa. As expected, PAM-based marine aerosols were associated with lower altitudes with respect to the sea level (96–325 m) as nearshore coastal areas are commonly characterized by lower elevations. MAR locations were concentrated along the coast of California, the Mediterranean Sea, and the Persian Gulf. MAR was rarely determined at island sites due to their typically low AOD values and the associated limitation for retrieving absorption-related aerosol properties (e.g., SSA). Aerosol sources within the BB aerosol type

(Figure 4) included sites typically influenced by wildfires (e.g., Bonanza Creek, Alta Floresta) or man-originated fires (e.g., Barcelona, Moscow and Skukuza). No BB models were present in Sudan or west Africa due to the lack of measurements, relatively high AOD, and/or major influence of dust. In Asia, BB sites were identified in Iran (e.g., IASBS), Russia (e.g., Yakutsk), Thailand, Burma, and Indonesia. Lastly, in Australia, BB aerosols were likely originated from wildfires of woodlands and forests (e.g., Lake Argyle and Lucinda).

AERONET sites for MIX were ubiquitous, and their location generally overlapped with the location of other aerosol types such as DU (e.g., Cairo_EMA), MAR (e.g., Arica), BB (e.g., Bonanza Creek), and U/I (e.g., UCLA). The AERONET sites classified as U/I corresponded with the largest cities around the world (e.g., Sao Paulo, Dragon Osaka Center, Mexico City). Also, PAM models labeled as U/I were present in highly populated areas in the U.S., Africa, Europe, and Asia (e.g., Fresno, Toulouse, Dalanzadgad, **Supplementary Appendix Table A3**).

In this study, an alternative evaluation of PAM's ability for classifying aerosol types in different geographic locations was based on MED comparisons. Independent AERONET measurements were used (i.e., remaining 25% of records not
A

<figure><img src="image_7.png" /></figure>

B

<figure><img src="image_8.png" /></figure>

C

<figure><img src="image_9.png" /></figure>

D

<figure><img src="image_10.png" /></figure>

**FIGURE 5** | Spectral variation of aerosol properties for selected PAM-derived clusters and associated types. Arithmetic average of normalized phase function at scattering angle θs = 0°) for forward panel (A) and backscattering angles [panel (B)], SSA [panel (C)], and g [panel (D)]. Curves with the same color symbol represent different aerosol models of the same type (Table 5). The numbers next to aerosol types represent PAM-derived subtypes (Table 5).

coincided with locations studied by Omar et al. (2005) in **Supplementary Appendix Table A3** (e.g., Venise, UCLA, Mexico City) and corresponded to sub-types 20 and 21. The identification of MIX-dominated sites based on PAM was not satisfactory as PAM predictions were only correct in two reference AERONET locations (Kanpur and Howland).

The AERONET sites labeled as DU by the PAM classifier coincided geographically with deserts (Supplementary Appendix Table A4 and Figure 4) (e.g., Mojave, Sonoran). To date, no AERONET sites exist in the Gobi Desert. However, its influence was manifested in western China. Predominance of dust along with observations made in Florida and the Caribbean Sea during summer was likely related to trans-Atlantic transport of dust originated in Sahara and Sub-Sahara (Sahel) deserts of Western Africa. As expected, PAM-based marine aerosols were associated with lower altitudes with respect to the sea level (96–325 m) as nearshore coastal areas are commonly characterized by lower elevations. MAR locations were concentrated along the coast of California, the Mediterranean Sea, and the Persian Gulf. MAR was rarely determined at island sites due to their typically low AOD values and the associated limitation for retrieving absorption-related aerosol properties (e.g., SSA). Aerosol sources within the BB aerosol type

(Figure 4) included sites typically influenced by wildfires (e.g., Bonanza Creek, Alta Floresta) or man-originated fires (e.g., Barcelona, Moscow and Skukuza). No BB models were present in Sudan or west Africa due to the lack of measurements, relatively high AOD, and/or major influence of dust. In Asia, BB sites were identified in Iran (e.g., IASBS), Russia (e.g., Yakutsk), Thailand, Burma, and Indonesia. Lastly, in Australia, BB aerosols were likely originated from wildfires of woodlands and forests (e.g., Lake Argyle and Lucinda).

AERONET sites for MIX were ubiquitous, and their location generally overlapped with the location of other aerosol types such as DU (e.g., Cairo_EMA), MAR (e.g., Arica), BB (e.g., Bonanza Creek), and U/I (e.g., UCLA). The AERONET sites classified as U/I corresponded with the largest cities around the world (e.g., Sao Paulo, Dragon Osaka Center, Mexico City). Also, PAM models labeled as U/I were present in highly populated areas in the U.S., Africa, Europe, and Asia (e.g., Fresno, Toulouse, Dalanzadgad, Supplementary Appendix Table A3).

In this study, an alternative evaluation of PAM's ability for classifying aerosol types in different geographic locations was based on MED comparisons. Independent AERONET measurements were used (i.e., remaining 25% of records not
<figure><img src="image_7.png" /><figcaption>**FIGURE 6** | Validation of aerosol models using MED in Solar Village. Predicted aerosol sub-types (A), MED<sup>min</sup> prediction performance (B). The mismatch between aerosol types is an agreement of 0. No blue circles when co-dominance is present.</figcaption></figure>

included in PAM analysis), and dominant clusters were
determined at the monthly scale for only four sites
(i.e., Solar Village, XiangHe, Mexico City, and Mongu).
Unfortunately, no MED<sup>min</sup> calculations were possible for
the MAR type due to the scarcity of measurements. In
general, the identification of aerosol types based on MED<sup>min</sup>
was very satisfactory for DU and BB case studies (section 2.4)
(i.e., Solar Village and Mongu, respectively). Indeed, the
annual mean of aerosol properties in PAM-derived DU-and
BB-related clusters allowed successful discrimination (i.e., 75
and 80% of months, respectively) of aerosol assemblages in
pre-selected reference sites. In Solar Village, DU-related sub-
types 4 and 9 were correctly predicted between January and
August, however, the MED<sup>min</sup> estimates between September
and December were poor due in part to a more diverse
composition of aerosols (e.g., co-dominance of U/I and
MIX) (Figure 6). The co-dominance among different DU
sub-types was also observed between January and March. In
Mongu, the dominant aerosol composition between June and
October was BB. The co-dominance of polluted continental
aerosols (i.e., U/I) during fall-winter (Omar et al., 2005) could
not be corroborated in this study due to the lack of
measurements (section 2.4). In XiangHe, there was no
correspondence between PAM and MED<sup>min</sup> classifications
of aerosol types as PAM-derived sub-types only matched in
9% of comparisons due to a strong co-dominance of aerosol

types different from MIX. However, MEDmin was more
sensitive as 67% of monthly data were labeled as MIX. The
validation results were also unclear in Mexico City, where
PAM clusters suggested a co-dominance between U/I and BB
and during different periods of the year. Also, MEDmin failed to
correctly identify the dominant aerosol composition as 90% of
monthly comparisons were associated with BB (e.g., sub-types
23 and 24). In this case, BB could be a dominant sub-type
during periods of moderate AOD (~0.4) for which the
AERONET retrieval is obtained. In contrast, our PAM
clustering is based on lower AOD values. As a result, the
validation procedure is dependent on the assumption
that the estimation of the aerosol type is conserved at
lower AOD.

Monthly variability of dominant clusters and their mean
aerosol properties are depicted in **Figure 7** for three AERONET
locations (GSFC, U.S. East coast; Kanpur, Northern India; and
Banizoumbou, West Africa). Most observations at GSFC were
made in spring-summer and showed a dominance of PAM
clusters linked to U/I and BB during spring-summer and fall,
respectively (**Figure 7A**). Winter months commonly exhibited a
more heterogeneous composition of aerosols with co-
dominance of one or more types (e.g., missing bars in
**Figure 7A**). The transition from spring-summer to fall was
associated with higher SSA(0.44) and lower AE<sub>e</sub> and FC<sub>v</sub> fine
values (**Figure 7B**). Also, the variability of AE<sub>a</sub> for both
dominant clusters was around 1. Months with missing mean
aerosol properties coincided with periods of the year where co-
dominance of aerosol sub-types occurred. The variation of
aerosol composition at Kanpur was characterized by MIX
and DU aerosol types during winter and spring-fall,
respectively (i.e., PAM clusters 15 and 2,5,8, respectively)
(**Figure 7C**). Unlike GSFC, the number of AERONET
observations was smaller and concentrated in the spring
months. Changes in aerosol properties of dominant clusters
coincided with an increase of AE<sub>a</sub> and SSA(0.44) from spring to
fall months(**Figure 7D**). Banizoumbou was characterized by
DU-related aerosol sub-types (e.g., 1,2,4 and 9) during all
seasons (**Figure 7E**). Likewise, DU-BB co-dominance was
mainly present between July and January. October was
associated with relatively high SSA(0.44) (~0.92), AE<sub>e</sub> (0.84)
and FC<sub>v</sub> fine (0.128) values (**Figure 7F**). Maximum mean values
of AE<sub>a</sub> (2.7) were calculated in June when AE<sub>e</sub> (0.1) and FC<sub>v</sub> fine
(0.08) were relatively low.

**Atmospheric Correction**

The performance of ACI and ACII for each AERONET-OC site is illustrated in **Table 6**. The performance metrics, and their ratios, suggest up to 30% improvements in estimating R<sub>rs</sub> (443) and R<sub>rs</sub> (482), corroborating the added value of PAM models. The differences in the performance of the two processors appear insignificant in the 561 and 655 nm bands, an observation consistent with the expected impact of absorbing aerosols (Nobileau and Antoine, 2005). The highest performance gain in R<sub>rs</sub> (443) is reported for the LISCO and Lake Erie sites. The LISCO site (~40 km northeast of New York City) is a region likely influenced by the U/I aerosols.
<figure><img src="image_11.png" /><figcaption>**FIGURE 6** | Validation of aerosol models using MED in Solar Village. Predicted aerosol sub-types (A), MED<sup>min</sup> prediction performance (B). The mismatch between aerosol types is an agreement of 0. No blue circles when co-dominance is present.</figcaption></figure>

included in PAM analysis), and dominant clusters were
determined at the monthly scale for only four sites
(i.e., Solar Village, XiangHe, Mexico City, and Mongu).
Unfortunately, no MED<sup>min</sup> calculations were possible for
the MAR type due to the scarcity of measurements. In
general, the identification of aerosol types based on MED<sup>min</sup>
was very satisfactory for DU and BB case studies (section 2.4)
(i.e., Solar Village and Mongu, respectively). Indeed, the
annual mean of aerosol properties in PAM-derived DU-and
BB-related clusters allowed successful discrimination (i.e., 75
and 80% of months, respectively) of aerosol assemblages in
pre-selected reference sites. In Solar Village, DU-related sub-
types 4 and 9 were correctly predicted between January and
August, however, the MED<sup>min</sup> estimates between September
and December were poor due in part to a more diverse
composition of aerosols (e.g., co-dominance of U/I and
MIX) (Figure 6). The co-dominance among different DU
sub-types was also observed between January and March. In
Mongu, the dominant aerosol composition between June and
October was BB. The co-dominance of polluted continental
aerosols (i.e., U/I) during fall-winter (Omar et al., 2005) could
not be corroborated in this study due to the lack of
measurements (section 2.4). In XiangHe, there was no
correspondence between PAM and MED<sup>min</sup> classifications
of aerosol types as PAM-derived sub-types only matched in
9% of comparisons due to a strong co-dominance of aerosol

types different from MIX. However, MEDmin was more
sensitive as 67% of monthly data were labeled as MIX. The
validation results were also unclear in Mexico City, where
PAM clusters suggested a co-dominance between U/I and BB
and during different periods of the year. Also, MEDmin failed to
correctly identify the dominant aerosol composition as 90% of
monthly comparisons were associated with BB (e.g., sub-types
23 and 24). In this case, BB could be a dominant sub-type
during periods of moderate AOD (~0.4) for which the
AERONET retrieval is obtained. In contrast, our PAM
clustering is based on lower AOD values. As a result, the
validation procedure is dependent on the assumption
that the estimation of the aerosol type is conserved at
lower AOD.

Monthly variability of dominant clusters and their mean aerosol properties are depicted in **Figure 7** for three AERONET locations (GSFC, U.S. East coast; Kanpur, Northern India; and Banizoumbou, West Africa). Most observations at GSFC were made in spring-summer and showed a dominance of PAM clusters linked to U/I and BB during spring-summer and fall, respectively (**Figure 7A**). Winter months commonly exhibited a more heterogeneous composition of aerosols with co-dominance of one or more types (e.g., missing bars in **Figure 7A**). The transition from spring-summer to fall was associated with higher SSA(0.44) and lower AE<sub>e</sub> and FC<sub>v</sub> fine values (**Figure 7B**). Also, the variability of AE<sub>a</sub> for both dominant clusters was around 1. Months with missing mean aerosol properties coincided with periods of the year where co-dominance of aerosol sub-types occurred. The variation of aerosol composition at Kanpur was characterized by MIX and DU aerosol types during winter and spring-fall, respectively (i.e., PAM clusters 15 and 2,5,8, respectively) (**Figure 7C**). Unlike GSFC, the number of AERONET observations was smaller and concentrated in the spring months. Changes in aerosol properties of dominant clusters coincided with an increase of AE<sub>a</sub> and SSA(0.44) from spring to fall months(**Figure 7D**). Banizoumbou was characterized by DU-related aerosol sub-types (e.g., 1,2,4 and 9) during all seasons (**Figure 7E**). Likewise, DU-BB co-dominance was mainly present between July and January. October was associated with relatively high SSA(0.44) (~0.92), AE<sub>e</sub> (0.84) and FC<sub>v</sub> fine (0.128) values (**Figure 7F**). Maximum mean values of AE<sub>a</sub> (2.7) were calculated in June when AE<sub>e</sub> (0.1) and FC<sub>v</sub> fine (0.08) were relatively low.

**Atmospheric Correction**

The performance of ACI and ACII for each AERONET-OC site is illustrated in **Table 6**. The performance metrics, and their ratios, suggest up to 30% improvements in estimating R<sub>rs</sub> (443) and R<sub>rs</sub> (482), corroborating the added value of PAM models. The differences in the performance of the two processors appear insignificant in the 561 and 655 nm bands, an observation consistent with the expected impact of absorbing aerosols (Nobileau and Antoine, 2005). The highest performance gain in R<sub>rs</sub> (443) is reported for the LISCO and Lake Erie sites. The LISCO site (~40 km northeast of New York City) is a region likely influenced by the U/I aerosols.
<figure><img src="image_8.png" /><figcaption>**FIGURE 7** | Site-specific monthly variability of PAM-derived aerosols. GSFC (A, B), Kanpur (C, D), Banizoumbou (E, F), dominant aerosol sub-type ids (left panels, solid circles), types (bold legend) and their mean properties (right panels), co-dominance of clusters (left panels, empty circles), $N_{clus}^{max}$ is the number of records for the dominant cluster.</figcaption></figure>

To further illustrate the improved spectral patterns in retrieved $R_{rs}$, six selected spectra from our matchup dataset are illustrated in **Figure 8**. The spectra reconstructed from ACII (dotted lines) better agree with the AERONET-OC spectra in the blue bands, indicating enhancements in $R_{rs}$ estimates with respect to those obtained from ACI. While the largest source of mismatch (e.g., biases) in retrieved $R_{rs}$ (from both ACI and ACII) can be remedied by modifying MDNs, the PAM models are shown to properly mend the shape of the spectra, providing a better fit to *in situ* measured spectra. The measured AOD(440) and retrieved SSA(440) (Level 1.5 and 2.0) are also annotated to inform the atmospheric conditions at the time of satellite overpass.

# DISCUSSION

## Clustering Approaches

This study's categorical analysis of aerosols was based on PAM, an unsupervised classification technique. There is no consensus in the literature regarding the preferential use of unsupervised over supervised methods for obtaining different aerosol classes. The main advantage of using a supervised algorithm is that clustering results are not abstract (i.e., aerosol classes have already been assigned to an aerosol type). Conversely, categories computed from unsupervised clustering must be
<figure><img src="image_12.png" /><figcaption>**FIGURE 7** | Site-specific monthly variability of PAM-derived aerosols. GSFC (A, B), Kanpur (C, D), Banizoumbou (E, F), dominant aerosol sub-type ids (left panels, solid circles), types (bold legend) and their mean properties (right panels), co-dominance of clusters (left panels, empty circles), $N_{clus}^{max}$ is the number of records for the dominant cluster.</figcaption></figure>

To further illustrate the improved spectral patterns in retrieved $R_{rs}$, six selected spectra from our matchup dataset are illustrated in Figure 8. The spectra reconstructed from ACII (dotted lines) better agree with the AERONET-OC spectra in the blue bands, indicating enhancements in $R_{rs}$ estimates with respect to those obtained from ACI. While the largest source of mismatch (e.g., biases) in retrieved $R_{rs}$ (from both ACI and ACII) can be remedied by modifying MDNs, the PAM models are shown to properly mend the shape of the spectra, providing a better fit to *in situ* measured spectra. The measured AOD(440) and retrieved SSA(440) (Level 1.5 and 2.0) are also annotated to inform the atmospheric conditions at the time of satellite overpass.

# DISCUSSION

## Clustering Approaches

This study's categorical analysis of aerosols was based on PAM, an unsupervised classification technique. There is no consensus in the literature regarding the preferential use of unsupervised over supervised methods for obtaining different aerosol classes. The main advantage of using a supervised algorithm is that clustering results are not abstract (i.e., aerosol classes have already been assigned to an aerosol type). Conversely, categories computed from unsupervised clustering must be
<table><caption>**TABLE 6** | Per-site matchup statistics (Eq. 5–7) for the performance of our experimental AC processors (ACI and ACII). The number of matchups is provided in the last column.</caption><thead><tr><th rowspan="2"></th><th colspan="4">ACI</th><th colspan="3">ACII</th><th rowspan="2">Metric Ratio</th><th rowspan="2">N</th></tr><tr><th>Band (nm)</th><th>MdSA (%)</th><th>MSA (%)</th><th>RMSLE</th><th>MdSA (%)</th><th>MSA (%)</th><th>RMSLE</th></tr></thead><tbody><tr><td rowspan="4">LISCO</td><td>443</td><td>32.6</td><td>28.5</td><td>0.12</td><td>30.8</td><td>22.8</td><td>0.10</td><td>1.06</td><td>1.25</td><td>1.20</td><td rowspan="4">7</td></tr><tr><td>482</td><td>13.8</td><td>13.7</td><td>0.07</td><td>11.7</td><td>12.8</td><td>0.06</td><td>1.18</td><td>1.07</td><td>1.25</td></tr><tr><td>560</td><td>21.8</td><td>27.0</td><td>0.11</td><td>18.4</td><td>23.0</td><td>0.10</td><td>1.18</td><td>1.17</td><td>1.14</td></tr><tr><td>655</td><td>30.1</td><td>33.2</td><td>0.13</td><td>26.9</td><td>26.4</td><td>0.11</td><td>1.12</td><td>1.26</td><td>1.20</td></tr><tr><td rowspan="4">Lake Erie</td><td>443</td><td>27.8</td><td>34.1</td><td>0.14</td><td>21.4</td><td>29.9</td><td>0.13</td><td>1.30</td><td>1.14</td><td>1.09</td><td rowspan="4">8</td></tr><tr><td>482</td><td>27.0</td><td>29.8</td><td>0.12</td><td>21.7</td><td>28.2</td><td>0.12</td><td>1.25</td><td>1.06</td><td>1.02</td></tr><tr><td>560</td><td>20.8</td><td>19.5</td><td>0.1</td><td>20.0</td><td>19.5</td><td>0.1</td><td>1.0</td><td>1.0</td><td>1.0</td></tr><tr><td>655</td><td>58.2</td><td>52.4</td><td>0.21</td><td>58.1</td><td>53.5</td><td>0.21</td><td>1.0</td><td>0.98</td><td>0.97</td></tr><tr><td rowspan="4">MVCO</td><td>443</td><td>24.0</td><td>27.3</td><td>0.1</td><td>19.0</td><td>25.3</td><td>0.1</td><td>1.3</td><td>1.1</td><td>1.0</td><td rowspan="4">24</td></tr><tr><td>482</td><td>22.1</td><td>25.7</td><td>0.12</td><td>19.3</td><td>24.1</td><td>0.12</td><td>1.14</td><td>1.07</td><td>1.05</td></tr><tr><td>560</td><td>20.0</td><td>26.9</td><td>0.1</td><td>20.0</td><td>27.7</td><td>0.1</td><td>1.0</td><td>1.0</td><td>1.0</td></tr><tr><td>655</td><td>32.0</td><td>31.0</td><td>0.14</td><td>28.5</td><td>31.2</td><td>0.14</td><td>1.12</td><td>0.99</td><td>0.99</td></tr><tr><td rowspan="4">Venise</td><td>443</td><td>16.6</td><td>17.6</td><td>0.1</td><td>15.2</td><td>15.9</td><td>0.1</td><td>1.1</td><td>1.1</td><td>1.1</td><td rowspan="4">94</td></tr><tr><td>482</td><td>15.8</td><td>16.7</td><td>0.08</td><td>15.3</td><td>15.5</td><td>0.07</td><td>1.03</td><td>1.07</td><td>1.08</td></tr><tr><td>560</td><td>12.2</td><td>14.7</td><td>0.1</td><td>11.7</td><td>15.4</td><td>0.1</td><td>1.0</td><td>1.0</td><td>1.0</td></tr><tr><td>655</td><td>20.1</td><td>25.7</td><td>0.13</td><td>21.3</td><td>27.3</td><td>0.13</td><td>0.95</td><td>0.94</td><td>0.97</td></tr><tr><td rowspan="4">Wave CIS</td><td>443</td><td>17.9</td><td>19.7</td><td>0.1</td><td>18.1</td><td>20.1</td><td>0.1</td><td>0.98</td><td>1.08</td><td>1.0</td><td rowspan="4">40</td></tr><tr><td>482</td><td>12.6</td><td>17.9</td><td>0.11</td><td>10.8</td><td>17.9</td><td>0.10</td><td>1.16</td><td>1.00</td><td>1.05</td></tr><tr><td>560</td><td>18.9</td><td>22.7</td><td>0.1</td><td>18.4</td><td>21.5</td><td>0.1</td><td>1.0</td><td>1.1</td><td>1.0</td></tr><tr><td>655</td><td>27.7</td><td>29.7</td><td>0.14</td><td>29.0</td><td>29.6</td><td>0.14</td><td>0.96</td><td>1.00</td><td>1.02</td></tr></tbody></table>

interpreted a *posteriori* using reference information (i.e., aerosol sources). Here, the maximum number of clusters was defined based on the Silhouette index, a metric successfully applied in several studies for discriminating aerosols (Omar et al., 2013; Äijälä et al., 2017; Szkop and Pietruczuk, 2019). Omar et al. (2005) computed a $K_{max}$ = 6 using an iterative criterion based on minimum cluster size. Unlike our findings, the Omar et al. (2005) estimates for $K_{max}$ were probably in the low range as their number of samples and sites was much smaller. To derive $K_{max}$ Giles et al. (2012) utilized multiple filters (e.g., 2-D scatterplots, size thresholds). They suggested up to 8 aerosol categories that might be underestimated due to our study's analysis of a larger number of AERONET sites. Russell et al. (2014) proposed a $K_{max}$ of seven that was established based on pre-defined aerosol sources, the Mahalanobis distance probabilities, and thresholds on aerosol properties. The smaller number of clusters in Russell et al. (2014) with respect to our study may be related to the spatially constrained analysis of AERONET sites and the strong correlation between some aerosol properties used for clustering. Overall, the finer clustering scheme adopted in this study was driven by the limited relevant information content in $ρ_t$, a physical constraint over natural waters where $πR_{rs}/ρ_t$ is generally < 10% (IOCCG, 2010). This challenge requires highly accurate quantification of aerosol contribution for robust retrievals of $R_{rs}$.

## Aerosol Types and Associated Properties

The spatial and temporal global distribution of aerosol types followed different patterns related to large-scale atmospheric

circulation (e.g., DU) or local conditions (e.g., U/I, BB, and MAR) (section 3.2). The prevalence of dust in the northern hemisphere has been linked to the air mass advection from deserts (e.g., Gobi and Taklamakan in north-central Asia, Saharan in Western Africa) (Tanré et al., 2003; Eck et al., 2005; Giles et al., 2011). The largest contribution of BB-related aerosols occurs during summer in northern high latitude locations, which tends to shift near and south of the equator in the northern hemisphere fall (Eck et al., 2003; Eck et al., 2009; Eck et al., 2013; Tilstone et al., 2020). Another finding was a major fraction of U/I during summer, perhaps due to humid, stagnant, and hot conditions allowing for secondary aerosol formation and hygroscopic growth due to humidification such as cloud processing (Eck et al., 2012; Eck et al., 2014).

PAM-derived aerosol properties in Table 5 showed unusually high $AE_e$ values for cluster 10 (MAR). Although $AE_e$ values up to 2.2 have been rarely measured (e.g., Lanai, Bermuda and Kaashidhoo (Smirnov et al., 2003a)), $AE_e$ values for maritime environments are typically below 1 when the aerosol optical depth at 500 nm is less than 0.15 (Smirnov et al., 2003a). The anomalous high mean value of $AE_e$ in our marine-dominated aerosol cluster is probably attributable to a relatively high contribution of fine-mode particles with a land origin (e.g., smoke) at Tremiti, Italy, and Chao-Jou, Taiwan. This hypothesis is supported by the relatively smaller $FC_v$ fine (~0.15) of maritime aerosols (Smirnov et al., 2003a; Sayer et al., 2012a). Mean values of $r_v$ and $σ_v$ for aerosol mode fractions reported in (Sayer et al., 2012a) did not explain the differences with our results as their fine and coarse particle distributions were characterized by a smaller radius ($r_n$ =
<table><caption>**TABLE 6** | Per-site matchup statistics (Eq. 5–7) for the performance of our experimental AC processors (ACI and ACII). The number of matchups is provided in the last column.</caption><thead><tr><th rowspan="2"></th><th colspan="4">ACI</th><th colspan="3">ACII</th><th rowspan="2">Metric Ratio</th><th rowspan="2">N</th></tr><tr><th>Band (nm)</th><th>MdSA (%)</th><th>MSA (%)</th><th>RMSLE</th><th>MdSA (%)</th><th>MSA (%)</th><th>RMSLE</th></tr></thead><tbody><tr><td rowspan="4">LISCO</td><td>443</td><td>32.6</td><td>28.5</td><td>0.12</td><td>30.8</td><td>22.8</td><td>0.10</td><td>1.06</td><td>1.25</td><td rowspan="4">7</td></tr><tr><td>482</td><td>13.8</td><td>13.7</td><td>0.07</td><td>11.7</td><td>12.8</td><td>0.06</td><td>1.18</td><td>1.07</td></tr><tr><td>560</td><td>21.8</td><td>27.0</td><td>0.11</td><td>18.4</td><td>23.0</td><td>0.10</td><td>1.18</td><td>1.17</td></tr><tr><td>655</td><td>30.1</td><td>33.2</td><td>0.13</td><td>26.9</td><td>26.4</td><td>0.11</td><td>1.12</td><td>1.26</td></tr><tr><td rowspan="4">Lake Erie</td><td>443</td><td>27.8</td><td>34.1</td><td>0.14</td><td>21.4</td><td>29.9</td><td>0.13</td><td>1.30</td><td>1.14</td><td rowspan="4">8</td></tr><tr><td>482</td><td>27.0</td><td>29.8</td><td>0.12</td><td>21.7</td><td>28.2</td><td>0.12</td><td>1.25</td><td>1.06</td></tr><tr><td>560</td><td>20.8</td><td>19.5</td><td>0.1</td><td>20.0</td><td>19.5</td><td>0.1</td><td>1.0</td><td>1.0</td></tr><tr><td>655</td><td>58.2</td><td>52.4</td><td>0.21</td><td>58.1</td><td>53.5</td><td>0.21</td><td>1.0</td><td>0.98</td></tr><tr><td rowspan="4">MVCO</td><td>443</td><td>24.0</td><td>27.3</td><td>0.1</td><td>19.0</td><td>25.3</td><td>0.1</td><td>1.3</td><td>1.1</td><td rowspan="4">24</td></tr><tr><td>482</td><td>22.1</td><td>25.7</td><td>0.12</td><td>19.3</td><td>24.1</td><td>0.12</td><td>1.14</td><td>1.07</td></tr><tr><td>560</td><td>20.0</td><td>26.9</td><td>0.1</td><td>20.0</td><td>27.7</td><td>0.1</td><td>1.0</td><td>1.0</td></tr><tr><td>655</td><td>32.0</td><td>31.0</td><td>0.14</td><td>28.5</td><td>31.2</td><td>0.14</td><td>1.12</td><td>0.99</td></tr><tr><td rowspan="4">Venise</td><td>443</td><td>16.6</td><td>17.6</td><td>0.1</td><td>15.2</td><td>15.9</td><td>0.1</td><td>1.1</td><td>1.1</td><td rowspan="4">94</td></tr><tr><td>482</td><td>15.8</td><td>16.7</td><td>0.08</td><td>15.3</td><td>15.5</td><td>0.07</td><td>1.03</td><td>1.07</td></tr><tr><td>560</td><td>12.2</td><td>14.7</td><td>0.1</td><td>11.7</td><td>15.4</td><td>0.1</td><td>1.0</td><td>1.0</td></tr><tr><td>655</td><td>20.1</td><td>25.7</td><td>0.13</td><td>21.3</td><td>27.3</td><td>0.13</td><td>0.95</td><td>0.94</td></tr><tr><td rowspan="4">Wave CIS</td><td>443</td><td>17.9</td><td>19.7</td><td>0.1</td><td>18.1</td><td>20.1</td><td>0.1</td><td>0.98</td><td>1.08</td><td rowspan="4">40</td></tr><tr><td>482</td><td>12.6</td><td>17.9</td><td>0.11</td><td>10.8</td><td>17.9</td><td>0.10</td><td>1.16</td><td>1.00</td></tr><tr><td>560</td><td>18.9</td><td>22.7</td><td>0.1</td><td>18.4</td><td>21.5</td><td>0.1</td><td>1.0</td><td>1.1</td></tr><tr><td>655</td><td>27.7</td><td>29.7</td><td>0.14</td><td>29.0</td><td>29.6</td><td>0.14</td><td>0.96</td><td>1.00</td></tr></tbody></table>

interpreted a *posteriori* using reference information (i.e., aerosol sources). Here, the maximum number of clusters was defined based on the Silhouette index, a metric successfully applied in several studies for discriminating aerosols (Omar et al., 2013; Äijälä et al., 2017; Szkop and Pietruczuk, 2019). Omar et al. (2005) computed a $K_{max}$ = 6 using an iterative criterion based on minimum cluster size. Unlike our findings, the Omar et al. (2005) estimates for $K_{max}$ were probably in the low range as their number of samples and sites was much smaller. To derive $K_{max}$ Giles et al. (2012) utilized multiple filters (e.g., 2-D scatterplots, size thresholds). They suggested up to 8 aerosol categories that might be underestimated due to our study's analysis of a larger number of AERONET sites. Russell et al. (2014) proposed a $K_{max}$ of seven that was established based on pre-defined aerosol sources, the Mahalanobis distance probabilities, and thresholds on aerosol properties. The smaller number of clusters in Russell et al. (2014) with respect to our study may be related to the spatially constrained analysis of AERONET sites and the strong correlation between some aerosol properties used for clustering. Overall, the finer clustering scheme adopted in this study was driven by the limited relevant information content in $ρ_t$, a physical constraint over natural waters where $πR_{rs}/ρ_t$ is generally < 10% (IOCCG, 2010). This challenge requires highly accurate quantification of aerosol contribution for robust retrievals of $R_{rs}$.

## Aerosol Types and Associated Properties

The spatial and temporal global distribution of aerosol types followed different patterns related to large-scale atmospheric

circulation (e.g., DU) or local conditions (e.g., U/I, BB, and MAR) (section 3.2). The prevalence of dust in the northern hemisphere has been linked to the air mass advection from deserts (e.g., Gobi and Taklamakan in north-central Asia, Saharan in Western Africa) (Tanré et al., 2003; Eck et al., 2005; Giles et al., 2011). The largest contribution of BB-related aerosols occurs during summer in northern high latitude locations, which tends to shift near and south of the equator in the northern hemisphere fall (Eck et al., 2003; Eck et al., 2009; Eck et al., 2013; Tilstone et al., 2020). Another finding was a major fraction of U/I during summer, perhaps due to humid, stagnant, and hot conditions allowing for secondary aerosol formation and hygroscopic growth due to humidification such as cloud processing (Eck et al., 2012; Eck et al., 2014).

PAM-derived aerosol properties in Table 5 showed unusually high $AE_e$ values for cluster 10 (MAR). Although $AE_e$ values up to 2.2 have been rarely measured (e.g., Lanai, Bermuda and Kaashidhoo (Smirnov et al., 2003a)), $AE_e$ values for maritime environments are typically below 1 when the aerosol optical depth at 500 nm is less than 0.15 (Smirnov et al., 2003a). The anomalous high mean value of $AE_e$ in our marine-dominated aerosol cluster is probably attributable to a relatively high contribution of fine-mode particles with a land origin (e.g., smoke) at Tremiti, Italy, and Chao-Jou, Taiwan. This hypothesis is supported by the relatively smaller $FC_v$ fine (~0.15) of maritime aerosols (Smirnov et al., 2003a; Sayer et al., 2012a). Mean values of $r_v$ and $σ_v$ for aerosol mode fractions reported in (Sayer et al., 2012a) did not explain the differences with our results as their fine and coarse particle distributions were characterized by a smaller radius ($r_n$ =
<figure><img src="image_9.png" /><figcaption>**FIGURE 8** | OLI-retrieved $R_{rs}$ shown for selected matchups. The spectra recovered through ACII (dotted) are more consistent with *in situ* (solid) than those estimated via ACI (dashed). The reported AOD(440) and SSA(440) (Level 1.5 and 2.0) are annotated for each instance. Note the differences in the ranges of y-axes.</figcaption></figure>

0.088 and 0.566 µm) and larger width (0.49 and 0.69) with respect to those associated to cluster 10.

The magnitude of $AE_e$ in cluster 14 was beyond the upper range of MIX aerosol assemblages, likely due to a greater contribution of particles derived from BB and U/I. In general, $AE_a$ values for MIX, U/I, MAR, and BB were low due to the increased dominance of black carbon particles. However, clusters 20 (U/I) and 25 (BB) had high $AE_a$ values indicating the possibility of brown carbon particles. Mean $AE_a$ values slightly less than one can be expected due to the uncertainty in calculating $AE_a$ from AERONET measurements (Giles et al., 2012).

## Spatial Comparisons

The spatial analysis of PAM-derived clusters revealed that aerosols of a continental origin mainly influence most AERONET sites (Figure 4). The marine aerosol type only represents ~ 4% of the entire database. This finding agrees with Omar et al. (2005), who suggested a minor marine contribution (i.e., 5%) based on nine years of measurements. Also, our results regarding BB contributions (26%) are comparable to those calculated by Omar et al. (2005) (26.6%).

Despite these similarities, our clustering results suggest that 1. natural sources (i.e., dust + marine) account for a more significant percentage of data due mainly to the underestimation of “dust” sites in Omar et al. (2005) (15%), 2. our U/I aerosol type contribution was 2-fold smaller, and 3. the contribution of the “mixed” aerosol type derived from PAM was more significant. Unlike our distribution of AERONET sites, the monitoring locations of Omar et al. (2005) did not include observations in Australia that were less dense in Asia, explaining the greater importance of “dust-type” particles in our study. The PAM classification generally only agreed with AERONET reference sites annually dominated by dust. Misclassification of DU was mainly due to the co-dominance of aerosol types during the same

month. This inconsistency was mainly associated with aerosol classes suggested by Omar et al. (2005). Also, inconsistencies in discriminating aerosols were attributed to dataset variations (e.g., sampling periods, aerosol load) at some AERONET sites chosen as reference aerosol sources (see Supplementary Appendix Table A3).

In our study, the Ilorin and Beijing sites corresponded to DU and U/I types (Supplementary Appendix Table A3), a result that apparently contradicts (Eck et al., 2010), who suggested a composition dominated by MIX. However, most AERONET records used in (Eck et al., 2010) were characterized by AOD(0.44) > 0.4 [see Figure 2 in Eck et al. (2010)] and different climatologies (e.g., only ten years for Ilorin). The accuracy of PAM for discriminating AERONET sites characterized by dominant aerosol types was also investigated in this study based on MED<sup>min</sup>, suggesting a satisfactory performance for detecting sites with a predominance of DU and BB types. Based on the k-means algorithm, Omar et al. (2005) also found that most comparisons associated with DU and BB were matched correctly. Unlike the annual dominant clustering approach used here, MED<sup>min</sup> had a better skill to separate MIX-dominated sites, however, this advantage was not observed for the U/I case study as aerosols were classified as BB. This interference was also detected by Omar et al. (2005) for three AERONET locations (Skukuza, IMS Metu Erdeml, and Dalanzadgad) and was attributed to smoke aging.

## Temporal Comparisons

Monthly changes on aerosol types at GSFC, Kanpur, and Banizoumbou were consistent with the aerosol type contributions computed in this study and local seasonal patterns of aerosol composition reported in the literature (Figure 7). The augmentation of smoke-related particles at GSFC during fall, dust events in Kanpur in spring, and seasonal alternate of DU-BB types in Banizoumbou are major
<figure><img src="image_13.png" /><figcaption>**FIGURE 8** | OLI-retrieved $R_{rs}$ shown for selected matchups. The spectra recovered through ACII (dotted) are more consistent with *in situ* (solid) than those estimated via ACI (dashed). The reported AOD(440) and SSA(440) (Level 1.5 and 2.0) are annotated for each instance. Note the differences in the ranges of y-axes.</figcaption></figure>

0.088 and 0.566 µm) and larger width (0.49 and 0.69) with respect to those associated to cluster 10.

The magnitude of $AE_e$ in cluster 14 was beyond the upper range of MIX aerosol assemblages, likely due to a greater contribution of particles derived from BB and U/I. In general, $AE_a$ values for MIX, U/I, MAR, and BB were low due to the increased dominance of black carbon particles. However, clusters 20 (U/I) and 25 (BB) had high $AE_a$ values indicating the possibility of brown carbon particles. Mean $AE_a$ values slightly less than one can be expected due to the uncertainty in calculating $AE_a$ from AERONET measurements (Giles et al., 2012).

## Spatial Comparisons

The spatial analysis of PAM-derived clusters revealed that aerosols of a continental origin mainly influence most AERONET sites (Figure 4). The marine aerosol type only represents ~ 4% of the entire database. This finding agrees with Omar et al. (2005), who suggested a minor marine contribution (i.e., 5%) based on nine years of measurements. Also, our results regarding BB contributions (26%) are comparable to those calculated by Omar et al. (2005) (26.6%).

Despite these similarities, our clustering results suggest that 1. natural sources (i.e., dust + marine) account for a more significant percentage of data due mainly to the underestimation of “dust” sites in Omar et al. (2005) (15%), 2. our U/I aerosol type contribution was 2-fold smaller, and 3. the contribution of the “mixed” aerosol type derived from PAM was more significant. Unlike our distribution of AERONET sites, the monitoring locations of Omar et al. (2005) did not include observations in Australia that were less dense in Asia, explaining the greater importance of “dust-type” particles in our study. The PAM classification generally only agreed with AERONET reference sites annually dominated by dust. Misclassification of DU was mainly due to the co-dominance of aerosol types during the same

month. This inconsistency was mainly associated with aerosol classes suggested by Omar et al. (2005). Also, inconsistencies in discriminating aerosols were attributed to dataset variations (e.g., sampling periods, aerosol load) at some AERONET sites chosen as reference aerosol sources (see Supplementary Appendix Table A3).

In our study, the Ilorin and Beijing sites corresponded to DU and U/I types (Supplementary Appendix Table A3), a result that apparently contradicts (Eck et al., 2010), who suggested a composition dominated by MIX. However, most AERONET records used in (Eck et al., 2010) were characterized by AOD(0.44) > 0.4 [see Figure 2 in Eck et al. (2010)] and different climatologies (e.g., only ten years for Ilorin). The accuracy of PAM for discriminating AERONET sites characterized by dominant aerosol types was also investigated in this study based on MED<sup>min</sup>, suggesting a satisfactory performance for detecting sites with a predominance of DU and BB types. Based on the k-means algorithm, Omar et al. (2005) also found that most comparisons associated with DU and BB were matched correctly. Unlike the annual dominant clustering approach used here, MED<sup>min</sup> had a better skill to separate MIX-dominated sites, however, this advantage was not observed for the U/I case study as aerosols were classified as BB. This interference was also detected by Omar et al. (2005) for three AERONET locations (Skukuza, IMS Metu Erdeml, and Dalanzadgad) and was attributed to smoke aging.

## Temporal Comparisons

Monthly changes on aerosol types at GSFC, Kanpur, and Banizoumbou were consistent with the aerosol type contributions computed in this study and local seasonal patterns of aerosol composition reported in the literature (Figure 7). The augmentation of smoke-related particles at GSFC during fall, dust events in Kanpur in spring, and seasonal alternate of DU-BB types in Banizoumbou are major
shifts in aerosol composition that agree with global temporal
patterns computed here (see Table 4). Omar et al. (2005)
concluded that aerosols at GSFC are the polluted continental
type with a higher contribution during summer as detected here.
However, Omar et al. (2005) pointed out that GSFC also
experiences periods of significant clean continental aerosol
(i.e., rural or MIX) from November to February, a finding not
supported by our classification that could be partly explained by
differences between datasets. Giles et al. (2011) showed that
aerosols at Kanpur have a complex composition with
absorption and spectral characteristics that resemble a MIX
aerosol type. Also, Giles et al. (2011) detected the
predominant influence of dust events during the spring
months. The year-round aerosol composition at Banizoumbou
was strongly influenced by dust originating from the Saharan
desert, as pointed out by Omar et al. (2005). Likewise, our models
support a switch in dominant aerosol composition occurring
between July-January (BB) and February-June (DU) (Rajot et al.,
2008).

Advantages and Disadvantages of Heritage
Models

PAM and A2010 models present major differences in terms of
aerosol properties (e.g., see Figures 1, 5). The range of r_n values
for the fine mode revealed that aerosol mixtures derived from the
clustering of AERONET measurements could have smaller
median diameters with respect to those suggested by A2010
models (Supplementary Appendix Table A5). Indeed, the
smallest r_n^fine values for the A2010 aerosol mixtures
(0.085 µm) were up to 2-fold larger than the minimum
median radius computed for PAM clusters. PAM-based
models have a smaller variability of r_n^coarse and intermediate
sizes (e.g., the radius was never above 1) with respect to A2010
models. Unlike the radius, the width of PSD was more variable for
fine and coarse fractions of PAM models as σ_n of A2010 models
was constant for each size-distribution mode (σ_n^fine = 0.437,
σ_n^coarse = 0.672). σ_n^fine (σ_n^coarse) for A2010 coincided with the
lower (upper) end of PAM-derived σ_n^fine (σ_n^coarse) values. The n_r
for A2010 models was more variable (1.379–1.528 and
1.371–1.518 for λ = 0.44 and 0.87 µm, respectively) with
respect to values associated with PAM clusters (1.417–1.502
and 1.425–1.1532). Despite this variation, the A2010 models
cannot capture the maximum values computed by PAM at λ =
0.87 µm. Unlike n_r, n_i values for A2010 models were less variable
and relatively small (e.g., 0.0022–0.010 and 0.00218–0.0099 for λ =
0.44 and 0.87 µm, respectively) with respect to those characterized
by PAM clusters (0.00132–0.052 and 0.00054–0.056). In terms of
aerosol types, two crucial questions arise:

1. What are operational aerosol models derived from A2010 not represented in nature (i.e., *in situ* information)?

2. How do microphysical properties of PAM-derived type MAR compare to those of RH-based models?

Regarding the first question, the following findings summarize
what aerosol characteristics are “missing” in A2010 models:

1. The median radius for the fine mode pointed out that ten DU, three MIX, and two U/I sub-types (i.e., clusters 1–10,13–15 and 18) are not part of A2010 models,

2. $\sigma_n^{\text{fine}}$ for the A2010 models only mimics the PAM cluster 19 with U/I type,

3. The range of $r_n^{\text{coarse}}$ values for A2010 can be constrained to a reduced number of theoretical models (2 model x 10 RH = 20 models) matching the PAM-derived range,

4. $\sigma_n^{\text{coarse}}$ for the A2010 models only resembles the PAM cluster 20 with U/I type,

5. The range of $n_r(0.44)$ values for A2010 can be constrained to 27 models matching the PAM-derived range,

6. Relatively high $n_r(0.87)$ values in PAM sub-types 2,6,14 and 20 are not reflected in A2010,

7. Relatively low and high $n_i(0.44)$ values (<0.002 and >0.01, respectively) in PAM models (eight U/I, four MIX, and three DU sub-types) are not presented in the A2010 models, and

8. Relatively low and high $n_i(0.87)$ values ($<0.0022$ and $> 0.099$, respectively) in PAM models (ten DU, six U/I, four MIX and four BB sub-types) are not included in any A2010 model.

Regarding the second question, the following can be
generalized (Supplementary Appendix Table A5):

1. The range of $r_n$, $\sigma_n$, and $n_r$ values used to compute A2010 models is narrower with respect to that characterizing PAM-derived MAR type.

2. The imaginary part of the refractive index is commonly underestimated by the A2010 aerosol models, and this difference is more remarkable at $\lambda = 0.44$ µm.

3. The arithmetic averages of $r_n^{\text{fine}}$, $r_n^{\text{coarse}}$, $n_r$, and $n_i$ for the PAM-derived MAR model represent A2010 models associated with an intermediate RH (i.e., RH = 75–80%).

4. The average width of fine mode (coarse) PSD was larger (smaller) for PAM-derived aerosol mixture.

Although the more realistic representation of σn and ni by PAM models, A2010 mixtures can be applied in a “dynamic” way, given the functionalities between RH and aerosol properties. Thus, the A2010 models may provide a finer temporal “tuning” as a function of humidity in coastal areas with marine-dominated aerosols. Also, RH-dependencies are still considered viable in marine-dominated locations not coinciding with AERONET sites.

Atmospheric Correction

The classification of aerosols in this study aimed to support an
enhanced representation of aerosol conditions for improved
global aquatic remote sensing science and applications. This
enhancement ultimately relies upon two important aspects: an
operational implementation of aerosol models and a realistic
representation of “natural” aerosol mixtures. NASA’s
operational AC method (Gordon and Wang, 1994), which
identifies aerosol types according to a NIR (or NIR-SWIR)
band-ratio, is generally ineffective when PAM-like aerosol
models are incorporated in the aerosol look-up-table
(Gordon, 1997). It is therefore imperative to devise
shifts in aerosol composition that agree with global temporal
patterns computed here (see Table 4). Omar et al. (2005)
concluded that aerosols at GSFC are the polluted continental
type with a higher contribution during summer as detected here.
However, Omar et al. (2005) pointed out that GSFC also
experiences periods of significant clean continental aerosol
(i.e., rural or MIX) from November to February, a finding not
supported by our classification that could be partly explained by
differences between datasets. Giles et al. (2011) showed that
aerosols at Kanpur have a complex composition with
absorption and spectral characteristics that resemble a MIX
aerosol type. Also, Giles et al. (2011) detected the
predominant influence of dust events during the spring
months. The year-round aerosol composition at Banizoumbou
was strongly influenced by dust originating from the Saharan
desert, as pointed out by Omar et al. (2005). Likewise, our models
support a switch in dominant aerosol composition occurring
between July-January (BB) and February-June (DU) (Rajot et al.,
2008).

Advantages and Disadvantages of Heritage
Models

PAM and A2010 models present major differences in terms of aerosol properties (e.g., see **Figures 1, 5**). The range of r<sub>n</sub> values for the fine mode revealed that aerosol mixtures derived from the clustering of AERONET measurements could have smaller median diameters with respect to those suggested by A2010 models (**Supplementary Appendix Table A5**). Indeed, the smallest r<sub>n</sub><sup>fine</sup> values for the A2010 aerosol mixtures (0.085 µm) were up to 2-fold larger than the minimum median radius computed for PAM clusters. PAM-based models have a smaller variability of r<sub>n</sub><sup>coarse</sup> and intermediate sizes (e.g., the radius was never above 1) with respect to A2010 models. Unlike the radius, the width of PSD was more variable for fine and coarse fractions of PAM models as σ<sub>n</sub> of A2010 models was constant for each size-distribution mode (σ<sub>n</sub><sup>fine</sup> = 0.437, σ<sub>n</sub><sup>coarse</sup> = 0.672). σ<sub>n</sub><sup>fine</sup> (σ<sub>n</sub><sup>coarse</sup>) for A2010 coincided with the lower (upper) end of PAM-derived σ<sub>n</sub><sup>fine</sup> (σ<sub>n</sub><sup>coarse</sup>) values. The n<sub>r</sub> for A2010 models was more variable (1.379–1.528 and 1.371–1.518 for λ = 0.44 and 0.87 µm, respectively) with respect to values associated with PAM clusters (1.417–1.502 and 1.425–1.1532). Despite this variation, the A2010 models cannot capture the maximum values computed by PAM at λ = 0.87 µm. Unlike n<sub>r</sub>, n<sub>i</sub> values for A2010 models were less variable and relatively small (e.g., 0.0022–0.010 and 0.00218–0.0099 for λ = 0.44 and 0.87 µm, respectively) with respect to those characterized by PAM clusters (0.00132–0.052 and 0.00054–0.056). In terms of aerosol types, two crucial questions arise:

1. What are operational aerosol models derived from A2010 not represented in nature (i.e., *in situ* information)?

2. How do microphysical properties of PAM-derived type MAR compare to those of RH-based models?

Regarding the first question, the following findings summarize
what aerosol characteristics are “missing” in A2010 models:

1. The median radius for the fine mode pointed out that ten DU, three MIX, and two U/I sub-types (i.e., clusters 1–10,13–15 and 18) are not part of A2010 models,

2. σ<sub>n</sub><sup>fine</sup> for the A2010 models only mimics the PAM cluster 19 with U/I type,

3. The range of r<sub>n</sub><sup>coarse</sup> values for A2010 can be constrained to a reduced number of theoretical models (2 model x 10 RH = 20 models) matching the PAM-derived range,

4. σ<sub>n</sub><sup>coarse</sup> for the A2010 models only resembles the PAM cluster 20 with U/I type,

5. The range of n<sub>r</sub>(0.44) values for A2010 can be constrained to 27 models matching the PAM-derived range,

6. Relatively high n<sub>r</sub>(0.87) values in PAM sub-types 2,6,14 and 20 are not reflected in A2010,

7. Relatively low and high n<sub>i</sub>(0.44) values (<0.002 and >0.01, respectively) in PAM models (eight U/I, four MIX, and three DU sub-types) are not presented in the A2010 models, and

8. Relatively low and high n<sub>i</sub>(0.87) values (<0.0022 and > 0.099, respectively) in PAM models (ten DU, six U/I, four MIX and four BB sub-types) are not included in any A2010 model.

Regarding the second question, the following can be
generalized (**Supplementary Appendix Table A5**):

1. The range of r<sub>n</sub>, σ<sub>n</sub>, and n<sub>r</sub> values used to compute A2010 models is narrower with respect to that characterizing PAM-derived MAR type.

2. The imaginary part of the refractive index is commonly underestimated by the A2010 aerosol models, and this difference is more remarkable at λ = 0.44 µm.

3. The arithmetic averages of r<sub>n</sub><sup>fine</sup>, r<sub>n</sub><sup>coarse</sup>, n<sub>r</sub>, and n<sub>i</sub> for the PAM-derived MAR model represent A2010 models associated with an intermediate RH (i.e., RH = 75–80%).

4. The average width of fine mode (coarse) PSD was larger (smaller) for PAM-derived aerosol mixture.

Although the more realistic representation of σ<sub>n</sub> and n<sub>i</sub> by PAM models, A2010 mixtures can be applied in a “dynamic” way, given the functionalities between RH and aerosol properties. Thus, the A2010 models may provide a finer temporal “tuning” as a function of humidity in coastal areas with marine-dominated aerosols. Also, RH-dependencies are still considered viable in marine-dominated locations not coinciding with AERONET sites.

**Atmospheric Correction**

The classification of aerosols in this study aimed to support an
enhanced representation of aerosol conditions for improved
global aquatic remote sensing science and applications. This
enhancement ultimately relies upon two important aspects: an
operational implementation of aerosol models and a realistic
representation of “natural” aerosol mixtures. NASA’s
operational AC method (Gordon and Wang, 1994), which
identifies aerosol types according to a NIR (or NIR-SWIR)
band-ratio, is generally ineffective when PAM-like aerosol
models are incorporated in the aerosol look-up-table
(Gordon, 1997). It is therefore imperative to devise
optimization techniques (e.g., polynomial fitting, ML) to
leverage a wide array of aerosol types (Gordon et al., 1997;
Chomko and Gordon, 1998). Compared to the performance of
other ML processors validated in the literature (Fan et al., 2021;
Pahlevan et al., 2021a; Brockmann et al., 2016), our
preliminary results (Table 6 and Figure 8) suggest that
MDNs, once trained extensively, can offer a significant
boost in the quality of retrieved $R_{rs}$. In this study, we
demonstrated that when PAM models are incorporated into
MDNs training (i.e., the model learns the spectral behavior of
absorbing aerosols in relation to $\rho_t$), the performance
improvements are noticeable for $\lambda < 500$ nm, particularly in
areas likely influenced by absorbing aerosols. These
performance gains (<30%) will have implications for the
retrieval of Chla and IOPs, such as phytoplankton
absorption spectra ($a_{ph}$) that are pivotal to the
characterization of phytoplankton properties (Bricaud and
Stramski, 1990; Hoepffner and Sathyendranath, 1991;
IOCCG, 2014; Uitz et al., 2015; Xi et al., 2015; Wang et al.,
2016) across challenging freshwaters and nearshore coastal
regions. Using paired *in situ* $R_{rs}$ and $a_{ph}$, (Pahlevan et al.,
2021b) showed that high-performing retrieval algorithms lead
to ~20% uncertainties in $a_{ph}$ ($\lambda < 500$ nm). To achieve this level
of uncertainty, satellite-derived $R_{rs}$ are expected to contain
minimal uncertainties (likely < 10% (Pahlevan et al., 2021a));
hence, any reduction in uncertainties in $R_{rs}$ opens
opportunities for robust assessments of phytoplankton
communities from space. The PAM models will also allow
for an increase in the number of useable (or valid) pixels
currently excluded due to high aerosol loading that leads to
frequent invalid $R_{rs}$ retrievals in the blue (IOCCG, 2010).
Enabling such capability, however, requires the development
of a viable processor that operates well under optically thick
aerosol conditions. We nevertheless note that improved
treatment of aerosol contribution to is only one key
component of reducing the total uncertainty budget. Proper
handling of skylight reflection (Gilerson et al., 2018) and
characterization of adjacency effects (Sanders et al., 2001)
are other major contributors that deserve further research.
A lack of knowledge surrounding the vertical distribution of
aerosols and uncertainties in ancillary data are among other
weak links in AC frameworks (Frouin et al., 2019; Mélin et al.,
2022).

We recall that future satellite instruments (Werdell et al., 2019; Ustin and Middleton, 2021), which are anticipated to offer measurements within the short blue bands and UV region, will benefit even more from the PAM models once combined with the A2010 models. Adding commercial satellite images to the publicly available datasets (Vanhellemont and Ruddick, 2018) will expand research and application opportunities, underscoring the need for more comprehensive aerosol models. Our models will be especially beneficial to emerging scientific applications, such as mapping and monitoring water quality indicators in urban environments (Miao et al., 2020) or in other regions typically influenced by dust (Dube et al., 2015) or wildfire events (Blake et al., 2020; Kramer et al., 2020). With the global rise of

temperature and extended periods of drought, these events
are becoming widespread, demanding AC processors capable
of tackling changing aerosol compositions and contents.
Consequently, subsequent efforts may focus on designing
novel methodology for effectively incorporating the PAM-
like models in existing/future AC processors.

CONCLUSION

The primary motivation of the present study was to advance NASA's existing aerosol models (A2010) for remote sensing applications in inland and nearshore coastal environments. An accurate representation of aerosol conditions is fundamental to obtaining reliable $R_{rs}$ products critical for approximating water constituents' concentration, composition, and optical properties. Following the compilation of daily historic AERONET data (1993-2019), an unsupervised clustering technique, the Partition around Medoids (PAM), was applied to SSA(0.44), $AE_e$, and $AE_a$ to determine unique aerosol models that differ from those of A2010. The comparison of aerosol microphysical properties derived from PAM and A2010 models demonstrated the added value of our empirically derived aerosol mixtures. A summary of our main findings include:

1. The validation of PAM-derived aerosol types agreed well with those reported at GSFC and Kanpur, two intensely studied AERONET sites where aerosol mixtures are typically influenced by dust and biomass burning smoke.

2. Identifying MIX and U/I types based on PAM sub-types was not always consistent with the sites suggested in the literature and the independent data analyzed here. This apparent disagreement is mainly attributed to our data selection criterion (AOD(0.44) ≤ 0.4), the lack of data in previously reported databases, and contrasting aerosol typification reported by several studies for same reference sites (e.g., Shirahama).

3. Improving aerosol representation over A2010 models is essential in areas where the aerosol composition is predominantly DU, U/I, and BB.

4. A practical assessment of the utility of PAM models in the context of atmospheric correction verified enhanced retrievals in Landsat-8/OLI's blue bands. The performance gain is most anticipated across regions influenced by absorbing aerosols and within the UV and blue spectral bands.

We assert that integrating PAM models in existing and/or
future AC processors will reduce uncertainties in
downstream products over lakes, rivers, reservoirs, and
nearshore environments from the past (Landsat-5), current
(Landsat-9), and future (e.g., Landsat Next) satellite
instruments. The PAM models are particularly expected to
advance the quality and useability of satellite products for
studying and monitoring waterbodies near aerosol sources
(e.g., coastal cities, urban waters, Sub-Saharan African lakes,
South American and Western North American water
resources).
optimization techniques (e.g., polynomial fitting, ML) to
leverage a wide array of aerosol types (Gordon et al., 1997;
Chomko and Gordon, 1998). Compared to the performance of
other ML processors validated in the literature (Fan et al., 2021;
Pahlevan et al., 2021a; Brockmann et al., 2016), our
preliminary results (Table 6 and Figure 8) suggest that
MDNs, once trained extensively, can offer a significant
boost in the quality of retrieved $R_{rs}$. In this study, we
demonstrated that when PAM models are incorporated into
MDNs training (i.e., the model learns the spectral behavior of
absorbing aerosols in relation to $\rho_t$), the performance
improvements are noticeable for $\lambda < 500$ nm, particularly in
areas likely influenced by absorbing aerosols. These
performance gains (<30%) will have implications for the
retrieval of Chla and IOPs, such as phytoplankton
absorption spectra ($a_{ph}$) that are pivotal to the
characterization of phytoplankton properties (Bricaud and
Stramski, 1990; Hoepffner and Sathyendranath, 1991;
IOCCG, 2014; Uitz et al., 2015; Xi et al., 2015; Wang et al.,
2016) across challenging freshwaters and nearshore coastal
regions. Using paired *in situ* $R_{rs}$ and $a_{ph}$, (Pahlevan et al.,
2021b) showed that high-performing retrieval algorithms lead
to ~20% uncertainties in $a_{ph}$ ($\lambda < 500$ nm). To achieve this level
of uncertainty, satellite-derived $R_{rs}$ are expected to contain
minimal uncertainties (likely < 10% (Pahlevan et al., 2021a));
hence, any reduction in uncertainties in $R_{rs}$ opens
opportunities for robust assessments of phytoplankton
communities from space. The PAM models will also allow
for an increase in the number of useable (or valid) pixels
currently excluded due to high aerosol loading that leads to
frequent invalid $R_{rs}$ retrievals in the blue (IOCCG, 2010).
Enabling such capability, however, requires the development
of a viable processor that operates well under optically thick
aerosol conditions. We nevertheless note that improved
treatment of aerosol contribution to is only one key
component of reducing the total uncertainty budget. Proper
handling of skylight reflection (Gilerson et al., 2018) and
characterization of adjacency effects (Sanders et al., 2001)
are other major contributors that deserve further research.
A lack of knowledge surrounding the vertical distribution of
aerosols and uncertainties in ancillary data are among other
weak links in AC frameworks (Frouin et al., 2019; Mélin et al.,
2022).

We recall that future satellite instruments (Werdell et al., 2019; Ustin and Middleton, 2021), which are anticipated to offer measurements within the short blue bands and UV region, will benefit even more from the PAM models once combined with the A2010 models. Adding commercial satellite images to the publicly available datasets (Vanhellemont and Ruddick, 2018) will expand research and application opportunities, underscoring the need for more comprehensive aerosol models. Our models will be especially beneficial to emerging scientific applications, such as mapping and monitoring water quality indicators in urban environments (Miao et al., 2020) or in other regions typically influenced by dust (Dube et al., 2015) or wildfire events (Blake et al., 2020; Kramer et al., 2020). With the global rise of

temperature and extended periods of drought, these events
are becoming widespread, demanding AC processors capable
of tackling changing aerosol compositions and contents.
Consequently, subsequent efforts may focus on designing
novel methodology for effectively incorporating the PAM-
like models in existing/future AC processors.

CONCLUSION

The primary motivation of the present study was to advance NASA's existing aerosol models (A2010) for remote sensing applications in inland and nearshore coastal environments. An accurate representation of aerosol conditions is fundamental to obtaining reliable $R_{rs}$ products critical for approximating water constituents' concentration, composition, and optical properties. Following the compilation of daily historic AERONET data (1993-2019), an unsupervised clustering technique, the Partition around Medoids (PAM), was applied to SSA(0.44), $AE_e$, and $AE_a$ to determine unique aerosol models that differ from those of A2010. The comparison of aerosol microphysical properties derived from PAM and A2010 models demonstrated the added value of our empirically derived aerosol mixtures. A summary of our main findings include:

1. The validation of PAM-derived aerosol types agreed well with those reported at GSFC and Kanpur, two intensely studied AERONET sites where aerosol mixtures are typically influenced by dust and biomass burning smoke.

2. Identifying MIX and U/I types based on PAM sub-types was not always consistent with the sites suggested in the literature and the independent data analyzed here. This apparent disagreement is mainly attributed to our data selection criterion (AOD(0.44) ≤ 0.4), the lack of data in previously reported databases, and contrasting aerosol typification reported by several studies for same reference sites (e.g., Shirahama).

3. Improving aerosol representation over A2010 models is essential in areas where the aerosol composition is predominantly DU, U/I, and BB.

4. A practical assessment of the utility of PAM models in the context of atmospheric correction verified enhanced retrievals in Landsat-8/OLI's blue bands. The performance gain is most anticipated across regions influenced by absorbing aerosols and within the UV and blue spectral bands.

We assert that integrating PAM models in existing and/or
future AC processors will reduce uncertainties in
downstream products over lakes, rivers, reservoirs, and
nearshore environments from the past (Landsat-5), current
(Landsat-9), and future (e.g., Landsat Next) satellite
instruments. The PAM models are particularly expected to
advance the quality and useability of satellite products for
studying and monitoring waterbodies near aerosol sources
(e.g., coastal cities, urban waters, Sub-Saharan African lakes,
South American and Western North American water
resources).
# DATA AVAILABILITY STATEMENT

The original contributions presented in the study are included in the article/Supplementary Material. The aerosol models detailed in Table 5 can be accessed via https://github.com/STREAM-RS.

# AUTHOR CONTRIBUTIONS

MM: Formal analysis, Investigation, Visualization; NP: Conceptualization, Funding acquisition, Project administration; MM, P-WZ, BS: Methodology; NP, DG, J-CR: Supervision; MM, BS, NP: Validation; MM, NP, DG: Roles/Writing - original draft; MM, NP, DG, J-CR, RL, PJW, AS: Writing—review and editing.

# FUNDING

This work was funded under NASA ROSES contract #80HQTR19C0015, Remote Sensing of Water Quality

# REFERENCES

Adrian, R., O'Reilly, C. M., Zagarese, H., Baines, S. B., Hessen, D. O., Keller, W., et al. (2009). Lakes as Sentinels of Climate Change. *Limnol. Oceanogr.* 54, 2283-2297. doi:10.4319/lo.2009.54.6_part_2.2283

Ahmad, Z., Franz, B. A., McClain, C. R., Kwiatkowska, E. J., Werdell, J., Shettle, E. P., et al. (2010). New Aerosol Models for the Retrieval of Aerosol Optical Thickness and Normalized Water-Leaving Radiances from the SeaWiFS and MODIS Sensors over Coastal Regions and Open Oceans. *Appl. Opt.* 49, 5545-5560. doi:10.1364/ao.49.005545

Äijälä, M., Heikkinen, L., Fröhlich, R., Canonaco, F., Prévôt, A. S., Junninen, H., et al. (2017). Resolving Anthropogenic Aerosol Pollution Types–Deconvolution and Exploratory Classification of Pollution Events. *Atmos. Chem. Phys.* 17, 3165-3197. doi:10.5194/acp-17-3165-2017

Alikas, K., Kangro, K., Randoja, R., Philipson, P., Asuküll, E., Pisek, J., et al. (2015). Satellite-based Products for Monitoring Optically Complex Inland Waters in Support of EU Water Framework Directive. *Int. J. Remote Sensing* 36, 4446-4468. doi:10.1080/01431161.2015.1083630

Balasubramanian, S. V., Pahlevan, N., Smith, B., Binding, C., Schalles, J., Loisel, H., et al. (2020). Robust Algorithm for Estimating Total Suspended Solids (TSS) in Inland and Nearshore Coastal Waters. *Remote Sensing Environ.* 246, 111768. doi:10.1016/j.rse.2020.111768

Barnes, B. B., Hu, C., Kovach, C., and Silverstein, R. N. (2015). Sediment Plumes Induced by the Port of Miami Dredging: Analysis and Interpretation Using Landsat and MODIS Data. *Remote Sensing Environ.* 170, 328-339. doi:10.1016/j.rse.2015.09.023

Barreto, A., Cuevas, E., Granados-Muñoz, M.-J., Alados-Arboledas, L., Romero, P. M., Gröbner, J., et al. (2016). The New Sun-Sky-Lunar Cimel CE318-T Multiband Photometer - a Comprehensive Performance Evaluation. *Atmos. Meas. Tech.* 9, 631-654. doi:10.5194/amt-9-631-2016

Bishop, C. M. (1994). *Mixture Density Networks*. Birmingham: Aston University Birmingham, 25.

Blake, D., Nyman, P., Nice, H., D'Souza, F. M. L., Kavazos, C. R. J., and Horwitz, P. (2020). Assessment of post-wildfire Erosion Risk and Effects on Water Quality in South-Western Australia. *Int. J. Wildland Fire* 29, 240-257. doi:10.1071/wf18123

Brajard, J., Jamet, C., Moulin, C., and Thiria, S. (2006). Use of a Neuro-Variational Inversion for Retrieving Oceanic and Atmospheric Constituents from Satellite Ocean Colour Sensor: Application to Absorbing Aerosols. *Neural Networks* 19, 178-185. doi:10.1016/j.neunet.2006.01.015

element, and the USGS Landsat Science Team Award #140G0118C0011.

# ACKNOWLEDGMENTS

We are thankful to Zia Ahmad for his initial involvement and feedback throughout the course of this study. We are also grateful to the Principal Investigators (PIs) of the AERONET-OC sites, whose datasets aided in evaluating the usefulness of our developed aerosol models for the atmospheric correction of Landsat-8/OLI images. The PIs include Giuseppe Zibordi, Sherwin Ladner, Hui Fang, Heidi Sosik, Tim Moore, Steve Ruberg, Menghua Wang, Alex Gilerson, and Samir Ahmed.

# SUPPLEMENTARY MATERIAL

The Supplementary Material for this article can be found online at: https://www.frontiersin.org/articles/10.3389/frsen.2022.860816/full#supplementary-material

Bricaud, A., and Stramski, D. (1990). Spectral Absorption Coefficients of Living Phytoplankton and Nonalgal Biogenous Matter: A Comparison between the Peru Upwelling Area and the Sargasso Sea. *Limnol. Oceanogr.* 35, 562-582. doi:10.4319/lo.1990.35.3.0562

Brockmann, C., Doerffer, R., Peters, M., Kerstin, S., Embacher, S., and Ruescas, A. (2016). Evolution of the C2RCC Neural Network for Sentinel 2 and 3 for the Retrieval of Ocean Colour Products in normal and Extreme Optically Complex Waters. *ESASP* 740, 54.

Cao, F., and Tzortziou, M. (2021). Capturing Dissolved Organic Carbon Dynamics with Landsat-8 and Sentinel-2 in Tidally Influenced Wetland- Estuarine Systems. *Sci. Total Environ.* 777, 145910. doi:10.1016/j.scitotenv. 2021.145910

Chomkō, R. M., and Gordon, H. R. (1998). Atmospheric Correction of Ocean Color Imagery: Use of the Junge Power-Law Aerosol Size Distribution with Variable Refractive index to Handle Aerosol Absorption. *Appl. Opt.* 37, 5560-5572. doi:10.1364/ao.37.005560

De Keukelaere, L., Sterckx, S., Adriaensen, S., Knaeps, E., Reusen, I., Giardino, C., et al. (2018). Atmospheric Correction of Landsat-8/OLI and Sentinel-2/MSI Data Using iCOR Algorithm: Validation for Coastal and Inland Waters. *Eur. J. Remote Sensing* 51, 525-542. doi:10.1080/22797254.2018.1457937

Derimian, Y., Karnieli, A., Kaufman, Y. J., Andreae, M. O., Andreae, T. W., Dubovik, O., et al. (2006). Dust and Pollution Aerosols over the Negev Desert, Israel: Properties, Transport, and Radiative Effect. *J. Geophys. Res. Atmospheres* 111. doi:10.1029/2005jd006549

Dey, S., Tripathi, S. N., Singh, R. P., and Holben, B. (2004). Influence of Dust Storms on the Aerosol Optical Properties over the Indo-Gangetic basin. *J. Geophys. Res. Atmospheres* 109, D20211. doi:10.1029/2004jd004924

Dube, T., Mutanga, O., Seutloali, K., Adelabu, S., and Shoko, C. (2015). Water Quality Monitoring in Sub-saharan African Lakes: a Review of Remote Sensing Applications. *Afr. J. Aquat. Sci.* 40, 1-7. doi:10.2989/16085914.2015.1014994

Dubovik, O., Holben, B., Eck, T. F., Smirnov, A., Kaufman, Y. J., King, M. D., et al. (2002a). Variability of Absorption and Optical Properties of Key Aerosol Types Observed in Worldwide Locations. *J. Atmos. Sci.* 59, 590-608. doi:10.1175/ 1520-0469(2002)059<0590:voaaop>2.0.co;2

Dubovik, O., Holben, B. N., Lapyonok, T., Sinyuk, A., Mishchenko, M. I., Yang, P., et al. (2002b). Non-spherical Aerosol Retrieval Method Employing Light Scattering by Spheroids. *Geophys. Res. Lett.* 29, 54-544. doi:10.1029/ 2001gl014506

Dubovik, O., and King, M. D. (2000). A Flexible Inversion Algorithm for Retrieval of Aerosol Optical Properties from Sun and Sky Radiance Measurements. *J. Geophys. Res.* 105, 20673-20696. doi:10.1029/2000jd900282