Self-supervised learning with multimodal remote sensed maps for seafloor visual class inference

2.1. Seafloor data acquisition

Data from multiple sources can be used for benthic habitat classification (see Figure 1). These can be broadly split as follows:

• Remote sensed mapping priors

Remote sensing mapping methods can seamlessly cover spatial extents of several kilo-hectares to kilometres squared. Examples include acoustic bathymetry from MBES providing depth information, and the backscatter intensity from SSS and MBES illustrating the hardness of substrate. Satellite images (SI) and aerial images (AI) taken by satellites and airborne cameras can also show the visual appearance of the seafloor in shallow water. The resolution of such data modes varies between tens of centimetres to tens of metres depending on acquisition conditions but is typically insufficient to directly determine substrate and habitat types without some separate form of observation.

Substrate and habitat classes can be directly identified through recovery of physical samples or high-resolution images. Sampling is an inherently destructive process and is also time-consuming where typically gathering tens to just over a hundred samples to characterise a region (Neettiyath et al., 2021; Usui et al., 2017). Seafloor visual imaging is a non-destructive process. But due to the strong attenuation of light in water, seafloor visual imaging requires the use of artificial lighting and cameras that are operated within a few metres above the seabed (Bodenmann et al., 2017). Even with high-altitude imaging setups, observational footprints are limited to <100 hectares per 24 h of observation (Thornton et al., 2021). However, the high resolution of imaging surveys (<1 cm) makes them suitable for identifying different types of seafloor habitats, substrata and species. Several prior works have demonstrated semi-supervised classification of seafloor imagery achieving high F1 accuracy scores (Ojala et al., 2002; Yamada et al., 2022).

A major challenge for large-scale seafloor habitat characterisation is the large gap between the extent covered by remote sensing and sampling or in situ camera observations. Remote sensing datasets cover vast regions at high-resolution, they can have measurements across many channels and the same geo-location can be described using multiple sensing modalities. Managing the resulting high-dimensionality is challenging for direct use of classifiers, and various approaches have been investigated to reduce dimensionality while retaining information that is useful for habitat interpretation. To efficiently identify boundaries between the different classes in a mapped region, automated interpretation methods often first project remote sensed data to a lower-dimensional feature space that aims to retain sufficient information while reducing information redundancy.

2.2. Feature extraction

Methods for feature extraction fall into two categories (Yamada et al., 2021b):

• Feature engineering, where features designed by humans are manually selected and combined

In feature engineering, Grey Level Co-occurrence Matrices (GLCMs) were first proposed in Haralick et al. (1973) as an effective and generalisable method for textural feature extraction. GLCMs combine image-derived features such as contrast and correlation to provide compact representations of grey-scale images. Zelada Leon et al. (2020) used 64 GLCM features and investigated the accuracy and repeatability of various methods for seafloor habitat classification. Fanlin et al. (2021) combined five GCLM features and angular response curves to represent acoustic images for seafloor classification. Scale-Invariant Feature Transform (SIFT) (Lowe, 1999) is an alternative to GLCMs that is extensively used to extract features from images. SIFT identifies points within an image where the intensity stands out from the local neighbourhood over different spatial scales. The distribution of these points forms a set of features that can be used to describe the image. Speeded-Up Robust Features (SURF) were developed to improve the SIFT by using optimised integral operations (Bay et al., 2006). Other forms of textural features include Local Binary Patterns (LBPs), which encode local relationships between the central pixel and its neighbouring pixels (Ojala et al., 2002). While engineered features have proven to be highly effective across a wide range of tasks, their performance may degrade under challenging conditions, such as in low-contrast or non-rigid scenes. Additionally, some methods, such as GLCMs, still require feature optimisation through combinations, which can be inconvenient to process multiple datasets.

Feature learning eliminates the need for manual optimisation and selection of predefined features. Instead, it leverages learning techniques to automatically extract features that efficiently describe a training dataset. Convolutional neural networks (CNNs) are widely used to flexibly learn representations of features directly from images, without the need for explicit feature engineering. The principle of deep CNNs has been analysed theoretically in Wiatowski and Bölcskei (2017) and many CNN architectures and training methods have been developed and applied to feature extraction and image classification tasks. AlexNet, proposed by Hinton et al. (2012), demonstrated that an increased number of layers in the network architecture enabled more complex feature representations, outperforming shallow networks in benchmarking studies. Later, deeper network architectures such as VGGnet (Simonyan and Zisserman, 2014) and Residual Network (ResNet) (He et al., 2016) were proposed and showed further improvements for classification tasks. A drawback of larger networks is the increased computational cost for network training, and CNNs such as MobileNet have attempted to reduce the number of parameters within the CNN architecture (Howard et al., 2017), for use in applications where computing resources are limited, for example, mobile robotic and sensing applications. Similar efforts have investigated the structure of CNNs to achieve efficient scaling of architectures (EfficientNet) for optimised performance under different computational constraints (Tan and Le, 2019). Recently, CNNs have been applied for automated interpretation of seafloor data. Mahmood et al. (2018) classified images of coral into nine classes, demonstrating the potential of CNNs to model habitats according to taxonomic class boundaries used in conservation ecology. However, a drawback of using CNNs is that traditional supervised learning approaches require large volumes of human-labelled training data to achieve accurate results. This is because the supervised classifiers simultaneously learn image descriptors (i.e., latent representation spaces) and delineate class boundaries in a single training process. The former requires hundreds to thousands of labelled instances for effective training, which is limiting since generating human-labelled data is time-consuming. Transfer-learning attempts to limit the number of domain-specific human-labelled training instances by first training CNNs on large, generic image datasets (Deng et al., 2009; Everingham et al., 2010; Lin et al., 2014), and then fine-tuning networks based on a smaller number domain-specific training data (Weiss et al., 2016). However, the performance of this approach is limited when the appearance of data in the target domain differs from the training repositories (Yamada et al., 2023).

Self-supervised learning has achieved state-of-the-art performance by separating the processes of feature learning and delineation of the learned feature space through semi-supervised classification (Chen et al., 2020a). Self-supervision trains CNNs to learn the features that best describe a dataset without the need for human-labelled datasets, which is an advantage because unlabelled datasets are more numerous and accessible. Instead, these approaches optimise CNNs parameters based on intrinsic structures that exist within unlabelled datasets (Chen et al., 2020a, 2020). Self-supervised learners attempt to map similar instances of data to the same region of a latent representation, or feature space, making them robust to small input or network weight perturbations (Samuli and Timo, 2017). Typically, artificially perturbed instances of the same sample are provided as pairs to the CNNs, where model weights are tuned so that the predictions (i.e., latent representation) are similar across the perturbations (Preciado-Grijalva et al., 2022). Self-supervised training methods have also been developed to pair different samples that are likely to have similar appearances based on certain rules. Yamada et al. (2022, 2021b) introduced a distance parameter when processing seafloor visual image datasets, where seafloor images that were taken close to each other were assumed to share similar properties since the footprint of an image is relatively small to the broader changes in seafloor substrates and habitats. Although this assumption will not always be satisfied, for example, near habitat transition areas, previous studies (Grant et al., 2024; Yamada et al., 2021b) showed that the approach is robust as long as the visual characteristics generally change over spatial scales larger than the distance parameter used.

Semi-supervised classifiers use intrinsic structures in the learned latent-representation space (Yamada et al., 2022, 2023) to reduce the number of labelled instances needed to delineate human class boundaries. In addition, this also reduces the susceptibility to human bias during training. For instance, Bijjahalli et al. (2023) developed a self-supervised learning framework based on variational autoencoders (VAE) to detect artificial objects in underwater imagery. The method applies clustering to the latent representation space and identifies samples that are distant from cluster centres. These images were considered anomalous and potentially showing artificial objects. Semi-supervision through human confirmation of artificial objects further improved learning performance. In order to avoid overfitting a small set of human-labelled data, several studies have investigated machine-guiding to identify samples for human-labelling based on the distribution of data in the latent representation space. For instance, Yamada et al. (2023) achieved 85% of the classification accuracy of an equivalent supervised CNN trained with 10,000 expert labels, using just 40 machine-guided labels. The approach used hierarchical K-means sampling to identify cluster representative images in the latent representation space, where algorithmically selected representative images improved the performance of classifiers compared to an equivalent number of randomly selected images. This approach was also found to improve performance under class imbalance by selecting the same number of images from each cluster found in the data.

2.3. Multimodal data interpretation

Combining different sensing modalities that capture diverse aspects of a measurement target has the potential to improve classification performance, and has been investigated in many applications, including seafloor habitat mapping (Rao et al., 2017), agriculture (Kang et al., 2023) and facial recognition (Nandi et al., 2022). Various approaches have been demonstrated, where these can be broadly grouped as early fusion, late fusion and middle fusion (Boulahia et al., 2021; Hong et al., 2020).

Early fusion merges raw or pre-processed sensor data before feature extraction. This ensures that joint features can be captured from multiple modalities at the onset of training, maximising the information derived from multimodal datasets (Cui et al., 2021; Gadzicki et al., 2020). A limitation of early fusion is an inherent sensitivity to spatial or temporal misalignment, which can arise from sensor resolution mismatches, geometric distortions, or changes in the environment between the acquisition of the different data modes. Boulahia et al. (2021) proposed early fusion to combine RGB images, depth maps and skeletal sequences for recognition of human actions. To address resolution mismatches between data modes, depth maps were resized to match RGB image sizes prior to feature extraction for downstream classification. Early fusion was also used to combine visual and thermal images for improved weed detection in agricultural applications (Zamani and Baleghi, 2023), where pre-processed visual images and thermal images were directly stacked ahead of feature extraction. In Jain et al. (2022), the authors demonstrated early fusion with contrastive self-supervised learning to interpret different satellite imaging modalities as similar pairs from the same geographic location. The results showed improved performance compared to a single mode when using multiple modalities. However, positional and geometric inconsistencies in satellite images are relatively limited compared to the remote sensing modes used in the marine domain, where different measurement physics and observation ranges are used, with inherent limitations in subsea positioning accuracy (Paull et al., 2014). Another factor is the longer time intervals between data acquisition due to mobilisation and survey logistics of marine operations, where even well-studied regions of the seafloor may only be visited once a year or less. The impact of these intermodal inconsistencies on multimodal learning is not understood.

Late fusion approaches combine the outputs from classifiers that are tailored to each sensing modality using a decision tree (Maki et al., 2011) or combine features derived from individual modalities by concatenating, or down-selecting a mixture of informative features as inputs for a classifier (Neettiyath et al., 2021; Takahashi et al., 2023). Advantages are that this approach is modular, where the introduction of new sensing modalities and corresponding mode-specific feature extractors does not impact other sensing modes, and there is flexibility as both feature engineering (Maki et al., 2011; Neettiyath et al., 2021) and feature learning (Takahashi et al., 2023) can be used. This approach can address issues of resolution mismatch that are common in remote sensing data as the feature extraction process can be coordinated between data modes. However, this also limits scalability as each data mode increases the computation cost and memory requirements. In addition, intermediate features that are potentially beneficial for classification may be missed, and conversely features from data modes that do not provide useful information can introduce noise into the final classification result. In Gunes and Piccardi (2005), the authors studied visual emotion recognition from facial expressions and body gestures, fusing the modalities at the decision level. Although the comparison experiments revealed that the early fusion method achieved higher classification accuracy, the late fusion made the feature extraction more flexible, where models for one modality can also be directly used in future tasks or recombined with other modes without re-training.

Middle fusion takes features that are individually extracted from each sensing modality and uses mathematical operations to derive a new set of features that express multimodal information. The distinction from late fusion is that the features used for classification are not a subset or combination of features extracted from each data mode, but instead features that individually contain information from multiple sensing modes. This approach maintains the flexibility of late fusion while allowing redundant information to be removed before classification. The method to fuse features from different modes is an actively researched area, where studies have demonstrated Bayesian optimisation (Ramachandram et al., 2017), genetic algorithms (Whitley et al., 1990), principal component analysis (Zelada Leon et al., 2020) and t-distributed Stochastic Neighbor Embedding (t-SNE) (Takahashi et al., 2023). In the study of Rao et al. (2017), features from visual images and bathymetry were extracted separately in a gated deep-learning architecture and then an intermediate shared layer was constructed to perform data fusion. The intermediate shared layer contains features of both the visual image dataset and bathymetry so that it is capable of predicting visual image and bathymetric features. In Takahashi et al. (2023), features derived from chemical signatures and holographic images of marine particles were combined using t-SNE. A comparison with late fusion, where the features derived from chemical signatures and imagery were concatenated, showed improved classification accuracy with t-SNE-based middle fusion, demonstrating the advantage of being able to recompute features based on multimodal inputs. Karpathy et al. (2014) compared early, late and middle fusion in the study of large-scale video classification, where a low-resolution context stream and a high-resolution fovea stream (i.e., images with non-uniform resolution, optimised based on some image region criteria) were fused in different ways. The result revealed that the middle fusion performed better than early and late fusion for this application (Feng et al., 2020). However, studies have found that the performance of these approaches is sensitive to data modalities, classification targets and network architectures (Feng et al., 2020), with no consensus being reached that any one approach consistently outperforms the others.

For benthic habitat classification, early fusion has the advantage that it can fully consider the information given by multiple modalities while remaining computationally efficient since only a single feature extractor and classifier are needed. However, early fusion is potentially more sensitive to data inconsistencies than late and middle fusion. Although methods such as resampling can be used to address the resolution mismatch, inconsistencies such as geometric distortion, positional offsets, and changes in the environment between the acquisition of different modes still exist. Although approaches to improve consistency across data layers exist, we argue that such spatial and temporal inconsistencies are inherent between seafloor mapping modalities, and so it is valuable to understand their impacts and develop robust approaches that minimise performance degradation.

2.4. Classification

In order to evaluate the feature learning performance, it is imperative to train a classifier to derive the ultimate predictive class labels. A significant challenge encountered during this process is data imbalance, a prevalent issue due to the heterogeneous distribution of habitats on the seafloor, resulting in imbalanced observations of AUVs. To address this, recent methodologies have been developed, focusing on rectifying data imbalance. These include strategies such as under-sampling and over-sampling (He and Garcia, 2009), along with algorithmic approaches like boosting (Singh and Purohit, 2015) and bagging (Hasib et al., 2020). Resampling techniques aim to equalise the representation of data classes by either augmenting or diminishing the number of samples. Conversely, boosting and bagging enhance training effectiveness through methods such as the utilisation of sub-datasets for repeated training and the aggregation of multiple weak classifiers to fortify the overall strength of the model. In this study, we deploy the SMOTE (synthetic minority over-sampling technique) (Fernández et al., 2018), a specific method in resampling techniques, to deal with the data imbalance problem owning to its easier implementation and low computation load compared to other strategies (Chawla et al., 2002).

Labelled visual image datasets are utilised to establish the ground truth through geo-registration processes. However, the footprint of these visual images, denoted as A , is approximately only 1 m², which is significantly smaller than the areas covered by patches of environmental priors, substantially exceeding A . Consequently, an individual patch in the training data can have multiple habitat classes in various proportions. This constitutes a probabilistic multi-class classification problem (Guo and Wang, 2015). Previous studies have applied both probabilistic and non-probabilistic classifiers to address this issue. For example, Qian et al. (2010) deployed a Support Vector Machine (SVM) to perform the multi-classes classification on human activities where several SVM classifiers were separately trained for each class, with final decisions made through a voting mechanism. Probabilistic classifiers, such as Gaussian Process Classifiers (GPC) (Huang, 2011), Gaussian Process Regression (GPR) (El-Mahallawy and Hashim, 2013) and Bayesian Neural Network (BNN) (Chaudhari and Tiwari, 2004), are also widely applied in real-scenarios because of their capability to give the possibilities of predictions. The probabilistic and non-probabilistic classifiers were extensively studied but the performances highly relied on the classification tasks (Chen et al., 2009). Therefore, both probabilistic and non-probabilistic classifiers are deployed in this study to thoroughly evaluate the self-supervised feature learning framework, aiming to eliminate the classifier-induced bias.

Deep learning techniques are also often used for direct classification. Both CNNs (Chaganti et al., 2020; Lee and Kwon, 2017) and transformers (Bhojanapalli et al., 2021) can be used to encode images into a feature space, where a MLP (multi-layer perceptron) can be added after these to predict class labels. Recently, transformers have demonstrated state-of-the-art performance in classification tasks (Han et al., 2022). Typically such approaches are trained directly through supervised learning. However, this requires large and well curated datasets (Zhou et al., 2021). Given the challenges of creating such a dataset for diverse marine habitats, supervised CNNs and transformer-based classifiers are not explored in this paper.

3.1. Convolutional early fusion

We leverage the proximity assumption to perform early data fusion using convolutional windows. Let I _s denote different data modes, for example, bathymetry and SSS. In our approach, the surveyed region is divided into geo-referenced patches of equal size, which we call y. The patch size determines the resolution of the habitat maps generated, and must be chosen to satisfy the following constraints:

• Contain multiple pixels to allow spatial patterns in the remote sensed data modes to be analysed

For the first point, although methods such as super-resolution can potentially be applied, we argue that this does not increase the basic information content. To fuse priors with different resolutions and channel depths, we first determine an appropriate patch size using the criteria below. The patch edge length r _d (in metres) must lie in the range:

r d ≥ σ pos

(1)

r d ≥ 10 ⌈ n min r res 10 ⌉

(2)

r d ≤ r hab

(3)

Each patch is discretised to match the feature learning network’s input layer, where we use zero-order sampling to match the resolution of each prior. This ensures the same geo-locations are indexed uniformly across each prior. A learnable 3 × 3 kernel, K _s , is passed over each patch to fuse the different priors, where for priors with multiple channels (e.g., satellite imagery with RGB layers), each channel is weighted equally. To populate each patch, data from each mode s is fused as follows (Gunes and Piccardi, 2005):

y ( N , E ) = ∑ s = 1 S ∑ m , n = − r d / 2 r d / 2 I s ( N + m , E + n ) ∗ K s ( m , n )

(4)

3.2. Self-supervised feature learning

3.2.1. Location-guided autoencoder (LGA)

The LGA is a self-supervised method that extends the autoencoder (AE) to consider location information during feature learning. The AE is a learning architecture that consists of two neural networks, an encoder and a decoder. The former maps the input, y, to a latent representation space h = f(y). The latter reconstructs the input from the latent representation as y _r = g(h). The AE aims to minimise the reconstruction error between y and y _r . The AE loss function is:

min ϕ , θ L rec = min ∑ y − y r 2

(5)

where ϕ and θ represent the weight and bias parameters of the encoder and decoder networks, respectively.

The LGA implements the proximity assumption by modifying the autoencoder loss function to consider location information (Yamada et al., 2021a). We first define p _ij , which relates the distance between two inputs y _i and y _j , for i ≠ j as:

p j ∣ i = exp − γ i − γ j 2 / 2 r loc 2 ∑ k ≠ i exp − γ i − γ k 2 / 2 r loc 2 ,

(6)

p i j = p j ∣ i + p i ∣ j 2 N b

(7)

where p _ij = 0 when i = j. γ is the location of the fused remote sensing patch. The distance parameter r _loc is a normalising factor for γ that defines the range over which location proximity, and therefore feature similarity is assumed. N _b is the minibatch size used for parameter optimisation in each training epoch.

The affinity q _ij is derived from h and is optimised based on p _ij . For q _ij when i ≠ j, it is defined by the Student’s t-distribution as:

q i j = 1 + h i − h j 2 − 1 ∑ k ≠ l 1 + h k − h l 2 − 1 ,

(8)

where q _ij = 0 for i = j, and h _i and h _j are features extracted from the early fused multimodal patches γ _i and γ _j . By defining the affinity matrices P and Q with p _ij and q _ij as their elements, the location loss L _loc is defined as the Kullback–Leibler (KL) divergence of P from Q:

L loc = K L ( P ‖ Q ) = ∑ i ≠ j p i j ⁡ log p i j q i j .

(9)

Minimising L _loc forces Q to approach P, which embeds the correlation between the early-fused data representations and the location metadata. The LGA loss function is modified as follows:

L = ( 1 − λ ) L rec + λ L loc

(10)

where λ weighs the location-based loss relative to the reconstruction loss. This can be considered a soft constraint since the t-distribution is heavy-tailed compared to Gaussian distributions, so when pairs of early fused data are initially far apart in the latent representation space, (i.e., their appearance is dissimilar), they are less strongly constrained by the regularisation term and can be flexibly embedded in different regions of the representation space. Since the loss function only loosely constrains autoencoder training based on probabilistic distributions, it is inherently robust to over-fitting location data. The distance over which similarity is assumed is a hyperparameter that can be tuned, where setting this to zero produces outputs that are identical to a standard autoencoder. Once the networks are trained, location information is no longer needed to extract features from the multimodal data.

The LGA in this study uses the AlexNet (Hinton et al., 2012) as the encoder architecture, and its inverted counterpart is used as the decoder. Convolutional layers are transformed into transconvolutional layers, and max pooling layers are converted to max unpooling layers. The workflow of the LGA is illustrated in Figure 3(a).OPEN IN VIEWER

Figure 3.

Self-supervised feature learning methods. Grey blocks indicate where parameters are optimised through the learning process. (a) The LGA introduces a geo-location-based loss term L _geo to the standard autoencoder. The geo-location loss term prioritises patterns found between multimodal patches that are geographically close to each other by minimising the KL divergence between location and feature affinity matrices. This imposes a soft location constraint that prioritises features that recur in nearby location. The distance over which similarity is assumed is a parameter that can be tuned, where setting this to zero produces outputs that are identical to a standard autoencoder. (b) GeoCLR modifies SimCLR by selecting positive (i.e., similar) pairs from nearby locations < r loc ), and a negative pairing from a random location ( > r loc ) . The contrastive loss function is used to learn feature embeddings from the multimodal inputs, which imposes a hard location constraint as inputs sampled from within the assumed similarity distance are forced to nearby regions of the feature space. When the similarity distance is set to zero, the positive pair is sampled from the same fused data input, and so the feature embeddings become identical to SimCLR.

3.2.2. GeoCLR

GeoCLR implements a hard location constraint through contrastive feature learning. It extends SimCLR, which learns feature representations by maximising agreement in the feature space between differently augmented views generated from the same input y (Chen et al., 2020a). In SimCLR, augmentations are randomly applied, where random crop followed by resizing back to the original size, random colour distortion and random Gaussian blur are commonly applied to generate positive pairs (y _i , y _j ), and a randomly selected input y _k , which is also randomly augmented and used as the negative pairing.

Features are extracted using a CNN encoder, h = f(y) and a projection head g(⋅) is used to map h to a smaller feature space, z = g(h), where the following contrastive loss function is applied:

ℓ i , j = − log exp sim z i , z j / τ ∑ k = 1 2 N b 1 [ k ≠ i ] ⁡ exp sim z i , z k / τ

(11)

where sim() represents cosine similarity, 1 [ k ≠ j ] ∈ { 0,1 } is the indicator function which is 1 if k ≠ i, and τ is the temperature parameter. The total loss of the minibatch samples is subsequently obtained as:

L = 1 2 N b ∑ k = 1 N b ℓ ( 2 k − 1,2 k ) + ℓ ( 2 k , 2 k − 1 )

(12)

where N _b is the minibatch size. In SimCLR (Chen et al., 2020a), the projection head is a small multi-layer perceptron with one hidden layer, and the feature encoding CNN is based on ResNet (He et al., 2016).

GeoCLR extends SimCLR and implements the proximity assumption by taking location information into account. Positive pairs are generated from distinct inputs [y _i , y _j ], where these are sampled from physically nearby locations subject to the following constraint:

( N i − N j ) 2 + ( E i − E j ) 2 ≤ r loc

(13)

where (N, E) are the northing and easting positions of the input data. The negative patch is selected from a random location outside the radius r _loc , and all the inputs are randomly augmented (Yamada et al., 2022). Once trained, location information is no longer needed to extract features from multimodal data. The GeoCLR workflow is shown in Figure 3(b). When the similarity distance is set to zero, the positive pair is sampled from the same fused data input, and so the feature embeddings become identical to SimCLR.

Since the contrastive loss embeds positive pairs to nearby regions of the feature space, GeoCLR can be considered as a hard location constraint. As r _loc increases, the positive pair y _i and y _j can be sampled from increasingly distant locations, which decreases the likelihood of them showing similar scenes. This can degrade the performance of the feature embedding, making GeoCLR inherently more sensitive to the choice of r _loc than the LGA. The GeoCLR in this study uses the ResNet18 CNN (He et al., 2016) as the encoder architecture.

3.3. Visual class sampling

We use semi-supervised learning to train classifiers to predict visual classes over the remote sensed region. The training and test labels used for benchmarking are generated from visual images (Massot-Campos et al., 2023; Yamada et al., 2022) or by human experts. Since patches generated to discretise the remote sensed priors are typically larger than the footprint of AUV gathered imagery, those that have overlapping imagery are likely to contain multiple visual class labels for training and testing, where these are not guaranteed to be from a single class (Figure 1). Previous studies have used the most frequent label (He and Garcia, 2009; Leevy et al., 2018), or the most central label (Rao et al., 2017) as a definitive single label to use for training and testing purposes. However, this neglects the mixing of different habitat classes near boundaries. To address this, we assign multiple labels to each patch, accompanied by their respective proportions determined based on their frequencies within each patch (Shields et al., 2020), which actually becomes a probabilistic multi-class classification problem.

3.4. Habitat classification and prediction

Classifiers are trained to correlate patterns between learned feature spaces and class labels derived from a random subset of the visual images at each site. Their performance is determined based on the F1 scores against reference class labels derived from the remaining visual images that were not used for classifier training. We note the use of location-based regularisation does not affect the validity of the results since location information is only used during training, and not during testing and inference.

We investigate four well-established probabilistic and non-probabilistic classifiers to delineate the visual class boundaries in the feature space. Compared to non-probabilistic classifiers, probabilistic classifiers provide uncertainty estimates that indicate confidence in their predictions. Regions with high uncertainty may benefit from additional observations to improve correlation between patch features and their corresponding visual habitat classes, but it is also possible that the remote sensing priors do not contain the information required to characterise the environment.

The Support Vector Machine (SVM) is chosen to represent a non-probabilistic classifier due to its proven robustness and performance over a variety of geospatial classification tasks (Ojala et al., 2002). For probabilistic classifiers, we compare the performance of the Gaussian Process Classifier (GPC), Gaussian Process Regression (GPR) and a Bayesian Neural Network (BNN). Based on a preliminary optimisation study using the grid search method (Syarif et al., 2016), we deploy the radial basis function kernel in the SVM, the rational quadratic kernel in GPC and the Matern kernel in GPR, respectively, as these were found to provide the best performance for each classifier. The BNN classifier uses a fully connected neural network using three layers with 64, 32 and 8 nodes, respectively. Classification accuracy scores are determined using the visual image-derived class proportions in each test patch for all classifiers.

4.1. Dataset description

Experiments are carried out on multimodal datasets from three marine protected areas (MPAs) in the UK: Darwin Mounds (DM), Studland Bay (SB) and Greater Haig Fras (HF), as shown in Figure 4. The MPAs have different substrate and habitat types, and the datasets consist of various remote sensed mapping data and seafloor imagery taken by different AUVs (see Table 1). Further information about the data modes used in our comparative experiments can be found in the supplemental material.OPEN IN VIEWER

OPEN IN VIEWER

Table 1.

Dataset Description. Multimodal Datasets Were Compiled From Three Different MPAs in the UK, Each With a Unique Habitat Type.

Site description
Name	Darwin mounds	Studland bay	Greater Haig Fras
Location	59.82 o N, 7.36 o W	50.65 o N, 1.94 o W	50.36 o N, 7.72 o W
Average depth	970 m	4 m	108 m
Habitat type	Cold-water coral	Seagrass	Mix rock, coarse and find sediments
Remote sensed mapping data, (resolution, m/pixel)
SSS	0.2	—	0.2
MBES backscatter	0.5	—	—
Bathymetry	2.0	2.0	4.0
Satellite images	—	3.0	—
Aerial images	—	0.34	—
AUV camera survey
Platform (camera system)	Autosub6000 (BioCam)	Smarty200	ALR(BioCam)
Imaging altitude, m	5	1	5
Resolution, mm	2.3	0.3	1.9
Trajectory, km	43	2.0	15.8
Observed area, h	15.3	0.2	7.9
Visual class labels (training: test)	12290:2252	4507:687	5831:1458
Number classes	4	7	4

4.2. Experimental setup

The patch size used in this study has r _d = 10 m, which is the smallest patch size that satisfies the criteria in equations (1)–(3) for our datasets. The patch window was shifted by an interval of 5 m and all patches were used for feature learning. Two sets of experiments were carried out to investigate the following aspects of performance:

• Effectiveness of multimodal feature learning with location regularisation for habitat classification

For the LGA and AE, AlexNet was used as the encoder network where the latent representation dimension was set to 16. The network was trained for 100 epochs with a learning rate of 1.0 × 10⁻⁵, λ = 1.0 × 10³ and the weight decay was set to 1.0 × 10⁻⁵, with a minibatch size of N _d = 128 following the recommendations of Yamada et al. (2023). For the LGA, the distance parameter was set to r _loc = 8 m to apply the soft proximity assumption to all neighbouring patches.

Experiments with GeoCLR and SimCLR used ResNet18 as the underlying architecture. The latent representation h and z dimensions were set to 512 and 16, respectively. The training epoch was set to 800, with a learning rate of 5.0 × 10⁻⁴, a weight decay of 1.0 × 10⁻⁴, temperature τ = 0.07 and minibatch size of N _d = 128, following the recommendation of Yamada et al. (2022). For GeoCLR the distance parameter was set to r _loc = 8 m to apply the hard proximity assumption to all neighbouring patches. The size of the latent representation space and distance parameter were set to the same value to ensure that results are comparable between the methods.

AE-LGA and SimCLR-GeoCLR network training epochs were fixed to 100 and 800, respectively. Network training took 7 h and 30 h, respectively, on a NVIDIA TITAN RTX 24 GB GPU. Computing the distance metrics from location metadata for LGA and GeoCLR took less than a minute, which is negligible compared to the total network training time.

Experiments were performed using both single mode and multimodal priors, where multimodal experiments combine remote sensed backscatter (SSS or MBES) or imagery (SI or AI), with bathymetric data (depth or relative depth (RD)). RD was calculated using equation (14):

R D ( N , E ) = D ( N , E ) − min D ( N , E ) d

(14)

To investigate the sensitivity of each method to the distance hyperparameter, we selected the remote sensed priors that give the best multimodal inference results and investigated the sensitivity of the results to tuning of the location regularisation (i.e., distance parameter r _loc ). The value of r _loc was varied between 0, 8, 15, 22, 44, 66, 100 m for both the LGA and GeoCLR, where r _loc = 0 corresponds to the AE and SimCLR, respectively. The final set of experiments investigates the robustness of the method to data inconsistencies between the remote sensed priors by introducing 2.5, 5, 7.5 and 10 m positional offsets between each mode being combined.

The macro F1 score, averaging precision for each predicted class and recall for each labelled class (Grandini et al., 2020), is used to evaluate multi-class prediction performance. Specifically, we compute the average class proportion error as the differences between the predicted and labelled class proportions for each pair of classes over the testing dataset. Based on this, a confusion matrix of accuracy can be constructed, enabling the calculation of the macro F1 score for probabilistic multi-class classification. For visualisation purposes, the habitat maps in Figures 7–9(e) and (f) assign the class with the highest predicted proportion per patch. Figures 7–9(g)–(j) show the proportion of an individual class within each patch.

4.3. Multimodal feature learning and classifier performance

Results comparing the performance of the different feature learning and classification methods are summarised in Figure 5, while full tables of results can be found in the supplemental material. The results show that early fused multimodal data consistently outperforms the use of a single mode, and that location regularisation improves the performance of both the AE based and contrastive feature learning models. The hard location constraint of GeoCLR outperforms the soft location constraint in LGA, with GeoCLR achieving 16.7% higher average F1 score overall conditions. The best results for each dataset are achieved by GeoCLR using multimodal data, reaching 79.0%, 56.8% and 83.0% for DM, SB and HF, respectively. The lower F1 score for the SB dataset is attributed to this being a harder classification problem, with a large number of visual classes that have a similar appearance (i.e., seagrass with varying percentage cover).OPEN IN VIEWER

For the DM dataset, combining SSS + RD achieves an F1 score of 79.0%, improving performance over the best-performing single mode SSS by 3.4%. Similar improvements are seen for SB, with AI + Depth achieving 56.8%, improving the best-performing single mode AI by 6.4%, and for HF SSS + Depth achieved 83.0%, outperforming the best single mode SSS by 6%. In general, GeoCLR improves F1 scores when combining modes, increasing the F1 score by an average of 5.27% over the best-performing individual mode it combines, demonstrating its ability to take advantage of multimodal patterns in early fused data. However, improvement is not guaranteed in cases where there is a large discrepancy in the information content of the layers being combined. For example, RD performs poorly for the SB and HF datasets (F1 scores <30% across all conditions), and when combined with higher-performing single modes, the combined F1 score falls short of the single mode score in some of the cases. However, even in these scenarios, the F1 score is significantly higher than the average of the modes being combined (on average 11.9% higher), indicating that rather than just naively merging information, feature learning actively rejects poor information content. This property is also demonstrated by LGA, which improves over the average of combined modes by 4.7%. An interesting results for the DM dataset is that although combining MBES + RD improves performance (65.6%) over each individual mode (58.5% and 62.0%), the overall performance is still significantly below SSS + RD (79.0%). The backscatter measurements of MBES originate from the same acoustic pulses used to determine RD (i.e., same sensor and signal, but using different processing). Therefore, these datasets have perfect alignment. On the other hand SSS backscatter data has an average positional offset of 12.0 m relative to RD due to localisation uncertainty and different measurement physics causing geometric distortions between the dataset (discussed in the supplemental material). Despite this the higher quality of acoustic backscatter information in the SSS outweighs any negative impacts of these spatial inconsistencies when combining the layers. A more detailed investigation into the effects of spatial inconsistencies will be discussed in the robustness analysis section.

Figure 6 summarises the performance of each classifier across all feature learners. Probabilistic classifiers (GPC, GPR, BNN) marginally outperform the non-probabilistic classifier (SVM) by 3.4% for equivalent feature and dataset conditions. The overall performance is less dependent on classifier choice than feature learners. However, probabilistic classifiers can provide uncertainty estimates for predictions, which may be advantageous in scenarios such as path planning where a larger proportion of observations could be gathered from such areas to improve the prediction certainty.OPEN IN VIEWER

4.4. Importance of location metadata

Figure 5 and Table 2 show LGA and GeoCLR improve performance compared to the equivalent AE and SimCLR, indicating a benefit from implementing the proximity assumption in feature learning. When using a single data mode, the average F1 scores using the LGA are 7.7% higher than AE, and 28.8% higher for GeoCLR than SimCLR under equivalent conditions, indicating significant improvements are achieved by incorporating location metadata, especially for the hard constraint implemented through contrastive learning. When using multimodal data, the relative improvements further increase to 8.8% and 37.8% for LGA and GeoCLR compared to AE and SimCLR, respectively. The latter result reflects that the AE and SimCLR feature learning are unable to take advantage of the multimodal data, where the AE only improves performance by 2.1% and SimCLR performance decreases by 2.0%. In contrast, LGA improves its performance by 3.2% and GeoCLR by 6.8% over the best-performing single mode. The effectiveness of LGA varies with dataset characteristics, exceeding AE F1 scores by 2.91%, 3.56% and 13.4% on the DM, SB and HF datasets on average, respectively. Notably, for the HF dataset the AE performance degrades with the use of multimodal inputs, whereas LGA improves significantly, demonstrating a fundamental difference in feature-learning behaviours. The results indicate that the proximity assumption and use of location metadata are critical to the enhancement of habitat classification performance when using multimodal data in feature learning. A more detailed assessment of the impact of location metadata on multimodal feature learning can be found in the supplemental material.

OPEN IN VIEWER

Table 2.

Average F1 Scores of Single Mode and Multi-Mode Priors Using AE, LGA, SimCLR and GeoCLR.

Dataset	Single mode				Multi modes
	AE	LGA	SimCLR	GeoCLR	AE	LGA	SimCLR	GeoCLR
DM	40.5 ± 1.6	43.7 ± 1.8	37.1 ± 2.1	59.7 ± 1.6	47.6 ± 1.4	49.6 ± 1.4	38.9 ± 2.0	67.6 ± 1.6
SB	18.7 ± 1.1	22.9 ± 2.6	16.7 ± 1.8	40.4 ± 2.3	23.3 ± 1.4	26.1 ± 2.2	21.8 ± 1.4	47.9 ± 3.0
HF	37.0 ± 1.3	42.2 ± 1.5	28.6 ± 2.0	51.9 ± 2.8	32.0 ± 1.4	57.8 ± 1.1	23.8 ± 1.8	75.9 ± 2.2

4.5. Habitat class prediction and class uncertainty

Figures 7–9 show the best performing single and multimodal remote sensing priors, visual class training inputs and classification results for the DM, SB and HF datasets, respectively. The predicted visual class distribution for single and multimodal scenarios corresponds to the best-performing feature learner and classifier pairing (see Table 3). In Figures 7 and 9 for the DM and HF datasets, panels (e) and (f) show the class with the highest predicted proportion within each patch, and panels (k) and (l) present uncertainty maps given by the classifier, which reflect the confidence of the classifier to predictions. Figure 8 also shows the same information for the SB dataset.OPEN IN VIEWER

OPEN IN VIEWER

Figure 8.

Habitat classification results for the SB dataset, showing (a) the figure legend; (b) the reference class distribution derived from visual images; (c) and (d) are the aerial images (AI) and Depth priors; (e) and (f) are classification maps inferred from AI and AI + Depth; (g) and (h) are example class proportion maps of class Rock/algae inferred from single mode and multimodal priors; (i) and (j) are example class proportion maps of Sediment inferred from single model and multimodal priors. (k) and (l) present the uncertainties given by the BNN classifier.

OPEN IN VIEWER

Figure 9.

Habitat classification results for the Greater Haig Fras, showing (a) the figure legend; (b) the used reference class distribution derived from visual images (full visual image dataset class distribution could be seen in supplemental materials); (c) and (d) are SSS and Depth data; (e) and (f) depict the habitat classification results derived from the single mode data SSS and multimodal data SSS + Depth, respectively; (g) and (h) are example class proportion maps of Bedrock inferred from single mode and multimodal priors; (i) and (j) are example class proportion maps of Sediment inferred from single mode and multimodal priors; (k) and (l) present the uncertainty in given by GPR classifiers.

OPEN IN VIEWER

Table 3.

Best-Performing Single and Multimodal Priors in Habitat Classification and the Corresponding Classifiers.

Dataset	DM	SB	GF
Single mode	SSS	AI	SSS
Multi-mode	SSS + RD	AI + Depth	SSS + Depth
Classifier	GPR	BNN	GPR

In the DM dataset, Figure 7(f), corresponding to multimodal inputs (SSS + RD), diagonal strips exist in classification maps, which are caused by the offsets seen between adjacent passes of the AUV when gathering bathymetry data. Although these are artefacts, the overall F1 achieved using multimodal data is higher than the best performing single mode input SSS (see Figure 7(e)), showing robustness to common practical issues relating to input data quality. The predicted proportion and distribution of Mound Top, where the majority of cold-water corals are found are relatively similar between the single and multi-mode scenarios in Figure 7(e) and (f), respectively.

The SB dataset results show a similar overall distribution pattern between Figure 8(f) for multimodal (AI + Depth) classification and Figure 8(e) for the best single mode (AI). The AI image has a boat visible in the data, where this region was masked during training and testing and adopted the most frequent predicted class of its nearest neighbours. Although the overall distribution patterns are similar, the Rock/algae in the bottom left predicted by the single mode Figure 8(e) is not predicted in the multimodal case in Figure 8(f), where the predictive uncertainty is also higher in this region compared to the single mode case. The multimodal individual class proportions show a transition from the Rock/algae class and the Seagrass 80%–100% class in this region.

For the HF dataset, we also see similar overall distribution patterns between the single mode (SSS) and multimodal (SSS + RD) inputs. However, Figure 9(l) shows less uncertainties in the Sediment class for the multimodal input compared to the single mode classification shown in Figure 9(k).

Overall, inconsistencies among multiple priors tend to increase predicted class and uncertainty. This is reflected by reduced F1 scores in these areas which is the expected behaviour and illustrates that the classifiers demonstrated in this work provide useful predictions of uncertainty in their class predictions. A detailed discussion of proportion maps for other classes is provided in the supplemental material.

4.6. Sensitivity to distance parameters

The parameter r _loc in both LGA and GeoCLR determines the distance over which the proximity assumption is applied during self-supervised feature learning. In order to understand its impact, a sensitivity analysis is conducted using the best priors and classifiers for each dataset (i.e., conditions in Table 3).

The results across the three datasets over an order of magnitude difference in the value of r _loc are shown in Figure 10. Compared to LGA, GeoCLR is more sensitive to the distance parameter, which is an expected outcome given the hard constraint implemented by contrastive learning. GeoCLR achieves the best results when location metadata is set so that only neighbouring patches are used as positive pairings, where the improvement deteriorates with the increasing r _loc , which can be explained as positive pairings being forced from regions of the data that are apart and so more likely to be of different visual classes. The results for single and multimodal inputs follow a similar trend with increasing r _loc from 8 to 100 m. However, even with over an order of magnitude of the r _loc value, GeoCLR outperforms the equivalent SimCLR that does not implement the location regularisation.OPEN IN VIEWER

The LGA is less sensitivity to variations in r _loc over the same range. This is expected since the soft constraint implemented via the modified loss function prioritises features of patches in geographically closer regions in a smooth manner, where the reconstruction loss term means that patches that do not have similar appearances in the first place can remain far apart in the feature space, minimising the risk of overfitting the proximity constraint. The results for single and multimodal inputs show a less obvious trend with increasing r _loc , where the performance gains over AE are generally smaller, and are maintained over the entire range of distance parameters tested.

4.7. Robustness analysis

Multimodal data has inherent positional offsets between individual data modes. Information about the actual positional offset magnitude and direction for multimodal data used in this study is given in the supplemental material (Table A1). Manually co-registered points between the data layers have large variability in their offset magnitude and direction, which can be explained as due to geometric distortions that are known to occur in some measurements modes, for examples, SSS. To investigate the robustness of our models to these practical aspects, we introduced artificial offsets in various directions (45°, 135°, 225° and 315°), with magnitudes of 2.5 m, 5 m, 7.5 m and 10 m, respectively, based on representative positional errors of 10 m in subsea mapping (Paull et al., 2014). No significant correlation was found between the directions of the true and artificial offset directions, and so here we present only bulk statistics for the different offset distance magnitudes.

The choice of patch size can influence the training and robustness of the model. Large patches can degrade performance, in particular around habitat transition areas where seafloor characteristics can change over small spatial scales causing nearby patches to have different characteristics. Conversely, smaller patches risk omitting significant patterns critical to classification. Therefore, determining an appropriate patch size requires careful consideration of prior resolution and the scale of habitats, as guided by equations (1)–(3).

The impact of introducing different artificial positional offsets on the performance of multimodal feature learning varies across the three datasets. As a reference, the results in Figure 11 include conditions that do not implement location regularisation during feature learning (i.e., r _loc = 0 corresponding to AE and SimCLR), and compare optimally tuned conditions r _loc = 8 m to r _loc = 100 m. Larger positional offsets gradually reduce the performance of multimodal inference but the overall reduction in performance remains < 5 % across all conditions. For all conditions where different positional offsets were introduced, LGA and GeoCLR continuously performed better than AE and SimCLR, which demonstrates the robustness of the method. As positional offsets increase from 2.5 m to 10 m, the performance of AE, LGA (r _loc = 8) and LGA (r _loc = 100) across three datasets decreases by an average of 5.0%, 3.8% and 3.4%, respectively, which shows that LGA leveraging location metadata could mitigate negative impacts induced by positional offsets, with slightly lower degradation seen when leveraging location regularisation over a large distance. However, the performance of SimCLR, GeoCLR (r _loc = 8) and GeoCLR (r _loc = 100) across three datasets decreases by 3.5%, 5.0% and 4.1%, where the location-guiding increases the sensitivity to positional offsets. However, we note that although the relative results are reduced, GeoCLR still achieved 30.4% higher F1 score than SimCLR, indicating overall benefits of considering location regularisation during feature learning.OPEN IN VIEWER

Based on the results, we can conclude that LGA is less sensitive to the distance parameter r_loc and more robust to positional offsets compared to GeoCLR, which implements a harder distance constraint in contrastive learning. However, it is important to note that GeoCLR remains the most effective feature learning method, outperforming AE, LGA and SimCLR in all the investigated scenarios.