Habitat Prediction of Bigeye Tuna Based on Multi-Feature Fusion of Heterogenous Remote-Sensing Data

: Accurate habitat prediction of Bigeye Tuna, the main fishing target of tuna pelagic fishery, is of great significance to the fishing operation. In response to the fact that most of the current studies use single-source data for habitat prediction, and the association between spatiotemporal features and habitat distribution is not fully explored and that this has limited the further improvement of prediction accuracy, this paper analyzes the spatiotemporal distribution of the characteristics of Bigeye Tuna’s highly migratory nature. Additionally, it puts forward a method of habitat prediction that utilizes heterosource remote-sensing data for the four-dimensional time–space–environment– spectrum (TSES) for deep-level feature extraction. First, a multi-source heterogeneous dataset was constructed by combining the spatiotemporal distribution characteristics of the product-level environmental remote-sensing data and the L1B-level original spectral remote-sensing data, and then a multi-branch, dynamic spatiotemporal feature extraction, Long Short-Term Memory Network (LSTM) time-series model was constructed to extract the characteristics of the heterogeneous data. This model was constructed to fully explore and utilize the multidimensional deep-level TSES distribution features affecting the habitat prediction. Finally, the two types of heterogeneous data were subjected to the weighted average-based decision-level fusion to obtain the final prediction results. The experimental results show that compared with other methods, the proposed method in this paper outperforms traditional machine-learning models and other single-source, data-based time-series models, with R 2 reaching 0.96278 and RMSE decreasing to 0.031361 in the validation experiments of these models. In contrast, the method in this paper demonstrates good generalization ability and achieves accurate prediction of future fishery distribution.


Introduction
The Bigeye Tuna is a large migratory fish belonging to the "Thunnidae" family, which is mainly distributed in the tropical and subtropical waters of the Atlantic, Indian, and Pacific Oceans.It is the most important species that is caught in longline tuna fisheries [1].The fishing and processing of the Bigeye Tuna provide an important export revenue for some coastal states.Due to its high market value, Bigeye Tuna fishery is a very promising pelagic fishery industry [2].It is important for the development and sustainability of Bigeye Tuna fisheries to know in advance the location of its central habitat and to plan fishing and production rationally because of Bigeye Tuna's highly migratory nature and the large migration distances that they travel during their different life cycles [3].
The Bigeye Tuna's habitat selection is a reflection of its migratory nature.Therefore, there is an inextricable relationship between the variability of the marine water environment in a given space and the spatiotemporal distribution of Bigeye Tuna habitats [4].At the same time, there is a distinct seasonal character to the migration of migratory fish.This characteristic is reflected in the distribution of habitats within future periods within adjacent seasons can be persistently affected by changes in the marine environment.Additionally, the distribution patterns of the habitats and their response to environmental factors within the same season exhibit notable similarities.The distribution of the Bigeye Tuna habitat shows a certain temporal pattern because of the seasonal characteristics.
In previous research on fishery prediction, the distribution of habitats in a given month was mainly predicted by studying the relationship between marine environmental factors and the formation of fishery habitats in that month.Spatial autoregression and spatial clustering methods were used by researchers such as Yuan et al. [5] to establish a habitat prediction model for Bigeye Tuna in the Indian Ocean and to explore the correlation between the distribution of the fishery and the sea surface height and between the sea surface temperature and the concentration of chlorophyll.The mean squared error of Yuan et al.'s model was around 0.1742.Meanwhile, a convolutional neural network CNN was proposed by Zhu et al. [6] to construct a prediction model for the Northwest Pacific soft-shelled habitat.This approach utilized tertiary sea surface temperature data obtained from MOIDS inversion to classify and predict high-and low-productivity fishing areas.Zhu et al.'s model's average accuracy could reach 66.6%.The effects of marine environmental factors on habitat formation in the time dimension are gradually being taken into account as research continues.For example, an integrated model based on a Random Forest, Bagging Decision Tree, C5.0 Decision Tree, Gradient Boosting Decision Tree, AdaBoost, and Stacking was established and proposed by Hou et al. [7] for the prediction of albacore tuna habitats in the South Pacific Ocean.Simultaneously, sea surface temperature and chlorophyll value data from one month before and after the study were included, resulting in a maximum prediction accuracy of 75.84%.The results of Hou et al.'s model were significantly improved compared to models that used only the current month's environmental data, which proves that environmental factors have a continuous effect on the distribution and growth of fish.In addition, the pattern of change in the distribution of fishery data itself is beginning to be utilized by some scholars for time-series prediction of small-scale fishing areas.The time-series analysis model ARIMA was applied by Song et al. [7] to predict the changes in the time series of small yellow croaker catches in the southern Yellow Sea.The results of Song et al.'s ARIMA model had an overall relative error of less than 14.74%, successfully predicted.Future habitat distribution based on yearly time-series changes in CPUE.This method is more practical than the traditional method of predicting the current month's habitat and is gradually becoming a new trend in the current fishery prediction.This analysis of previous research work reveals that most of the previous studies are limited to product-level marine remote-sensing data or sequence network models with single-source data as the input for prediction [8][9][10][11][12][13][14][15][16].However, the more homogeneous remote-sensing information of the marine environment is not conducive to realizing accurate analysis and prediction of habitat prediction of changes in the fishery.The effective utilization of original ocean remote-sensing data and the realization of accurate habitat prediction through heterogeneous data fusion have become possible with the continuous development of satellite remote-sensing technology and geographic information systems using high spatiotemporal resolution.A decision-level fusion model for feature extraction of heterogeneous data was proposed by Han et al. [17] in order to solve the problem of difficult fusion of heterogeneous sources and heterogeneous remote-sensing data features extraction.In this model, data on environmental factors and multi-spectrum data were extracted by a convolutional neural network, respectively, temporal features of spatiotemporal factors were extracted using LSTM structure, and finally, the decision-level fusion method was constructed to realize the effective fusion of multi-features of heterogeneous data.Although the method has initially demonstrated the feasibility of using remote-sensing data from heterogeneous sources to predict habitat distribution, an in-depth study of the factors influencing the distribution of monsoon migratory fish habitats is still lacking.Additionally, the method is currently limited to predicting the distribution of habitats for the current month and does not effectively validate the accuracy of future habitat predictions.
Therefore, in order to enhance the accuracy and effectiveness of Bigeye Tuna habitat prediction, this paper builds on the above research to propose an innovative decision-level fusion model for multi-feature extraction and integration based on different spatiotemporal dimensional features of heterosource remote-sensing data, constructing multi-source heterosource datasets that leverage the spatiotemporal distribution characteristics of environmental remote-sensing data and original optical remote-sensing data to fully explore the temporal, spatial, environmental, and spectral (Time-Space-Environment-Spectrum, TSES) four-dimensional feature information influencing habitat prediction and fishing grounds.This method constructs multi-branch time-series feature extraction models for the structural characteristics of different datasets and performs in-depth extraction and fusion of the four-dimensional features of TSES.In the product-level dataset of marine environmental factors, the correlation between marine environment information and CPUE, as well as the change patterns of CPUE, are explored in both yearly and monthly dimensions.Simultaneously, the ConvLSTM2D time-series model is constructed for in-depth extraction of marine environment features in two-dimensional space under fixed time dimension, and the spatiotemporal factor information of the habitat is extracted for feature fusion using the 1D-CNN model.In the datasets of multi-source remote sensing, the dynamic changes of different band information over time are mined using the LSTM model to map with CPUE in multi-source remote-sensing datasets, which is then combined with spatiotemporal factor information for heterogeneous data feature fusion.Finally, accurate prediction of Bigeye Tuna habitats is achieved through decision-level fusion based on weighted averaging.The experimental results show that the method proposed in this paper provides better predictions compared to using only product-level remote-sensing data or L1B-level original remote-sensing data.The accuracy of habitat prediction was further improved by the effective fusion of the four-dimensional feature information of TSES, enabling a deeper analysis of the patterns in habitat distribution and the spatiotemporal marine environment of Bigeye Tuna.The practical value of the habitat prediction was further enhanced by the temporal prediction of the fishery based on the three temporal information of fishery, environment, and spectrum, providing accurate forecasts of center fishery in the future.

Data Sources
The fisheries' data were obtained from the Indian Ocean Tuna Commission (IOTC, https://iotc.org/,accessed on 30 March 2022) longline statistics.A total of 3975 records of Bigeye Tuna in the Southwest Indian Ocean, covering the months of January and May to December from 2003 to 2013, were selected.The data included year, month, longitude, latitude, number of hooks, and catch mass (kg), with a temporal resolution of months and spatial resolution of 5 • × 5 • .The study area ranged from 2 The growth activities of Bigeye Tuna are influenced differently by various depths and types of marine environmental factors.These factors at specific depth ranges were selected as experimental data according to the habitat preference of Bigeye Tuna and the impact of different depths of marine environmental factors on their habitat selection [18,19], including five types of environmental data: 0~300 m sea surface temperature (SST), chlorophyll concentration (Chlα), sea level anomaly (SLA), 0~400 m sea surface salinity (SSS) and 200~400 m dissolved oxygen concentration (O).The environmental data were all obtained from the Copernicus Marine Service (https://marine.copernicus.eu/,accessed on 30 March 2022), spanning the months of January and May to December from 2003 to 2013, with a temporal resolution of months and a spatial resolution of 0.25 The marine environmental factor data are tertiary product-level data, which may lose some of the effective information between spectral channels in the process of generation, making it difficult to accurately determine the correlation between environmental factors and habitat location.To address this, the original L1B-class multispectral remote-sensing data from the MODIS-Aqua satellite, of which 8~16 bands and 31~32 bands are related to the marine environment, are added as the experimental data of this experiment after data preprocessing.These data span from January 2003 to December 2013, with the temporal resolution of days and the spatial resolution of 5 × 5 km.The original remotesensing data were obtained from the MODIS data download address provided by NASA (https://ladsweb.modaps.eosdis.nasa.gov/,accessed on 30 March 2022).

Fisheries' Data Processing
(1) Calculate CPUE Catch Per Unit Effort (CPUE) is a measure commonly used to reflect the abundance of a fishery resource.In this experiment, CPUE is calculated using the number of hooks and the mass of the catch in kilograms within a 5 • × 5 • grid (in kg/thousand hooks).The formulae are as follows: where CPUE (i,j) , F (i,j) , and H (i,j) denote the CPUE, monthly total catch quality, and monthly total number of hooks released for the range of longitude i and latitude j, respectively.
(2) Expansion of monthly data grids The original fisheries data are organized on a 5 • × 5 • grid, which is too coarse to capture the distribution of habitats in detail for each month.At the same time, fisheries' data for some months are partially missing, resulting in incomplete coverage of the whole study area and hindering the analysis of the distribution pattern of fish stocks in different months.The inverse-distance weighting (IDW) method was employed in this experiment to address the above issues and align the spatial resolution of the environmental data by interpolating and expanding the fisheries data on a monthly basis over the geographic extent of the studied habitats.
The inverse-distance weighting (IDW) method is one of the commonly used interpolation techniques in Geostatistics for regridding irregular points or gridded data over a spatial extent.The basic reasoning of IDW is that all known measurement points exert a local influence on the complementary points, with the influence inversely proportional to the distance between them-the closer the points, the greater the influence, and vice versa.If the number of known measurement points in the spatial extent that influence the complementary points is n, and the weights are assigned according to the distance between the measurement points and the complementary points, decreasing as the distance increases [20].The complementary points F are calculated as follows: where F denotes the result obtained by interpolation of the inverse-distance weight method, n denotes the number of measurement points that influence the interpolated point, f i denotes the value of the i-th measurement point, and ω i denotes the weight of the i-th measurement point.The formula for ω i is as follows: where p is the weight coefficient, which can be any positive integer, and the default value of 2 is chosen for this experiment; h i denotes the Euclidean distance between the i-th measurement point and the complementary point.
In this experiment, the fishery data were divided by month to study the changing pattern of the monthly habitat better, and the data of each month were regridded using the inverse-distance weighting method.The single-month habitat data for Bigeye Tuna were generated by dividing the area into 1 • × 1 • grid.The weights were assigned based on the distances from the grid points to the known points, and the full values of the filled grid points were calculated.For example, the grid expansion before and after regridding for January 2006 is shown in Figures 1 and 2.
where p is the weight coefficient, which can be any positive integer, and the default value of 2 is chosen for this experiment; hi denotes the Euclidean distance between the i-th measurement point and the complementary point.
In this experiment, the fishery data were divided by month to study the changing pattern of the monthly habitat better, and the data of each month were regridded using the inverse-distance weighting method.The single-month habitat data for Bigeye Tuna were generated by dividing the area into 1° × 1° grid.The weights were assigned based on the distances from the grid points to the known points, and the full values of the filled grid points were calculated.For example, the grid expansion before and after regridding for January 2006 is shown in Figures 1 and 2.  In this experiment, the gridded fisheries data were expanded to 329 entries per month by the inverse-distance weighting method, totaling 32,571 entries for January and May to December 2003-2013.The spatial resolution of the fisheries data was improved from 5° × 5° to 1° × 1°, providing a more detailed portrayal of the monthly habitat distribution of Bigeye Tuna.The increase in the number of samples reduced the likelihood of random errors in model training, while the gap in spatial resolution with product-level environmental data was narrowed, facilitating subsequent data-matching.where p is the weight coefficient, which can be any positive integer, and the default value of 2 is chosen for this experiment; hi denotes the Euclidean distance between the i-th measurement point and the complementary point.
In this experiment, the fishery data were divided by month to study the changing pattern of the monthly habitat better, and the data of each month were regridded using the inverse-distance weighting method.The single-month habitat data for Bigeye Tuna were generated by dividing the area into 1° × 1° grid.The weights were assigned based on the distances from the grid points to the known points, and the full values of the filled grid points were calculated.For example, the grid expansion before and after regridding for January 2006 is shown in Figures 1 and 2.  In this experiment, the gridded fisheries data were expanded to 329 entries per month by the inverse-distance weighting method, totaling 32,571 entries for January and May to December 2003-2013.The spatial resolution of the fisheries data was improved from 5° × 5° to 1° × 1°, providing a more detailed portrayal of the monthly habitat distribution of Bigeye Tuna.The increase in the number of samples reduced the likelihood of random errors in model training, while the gap in spatial resolution with product-level environmental data was narrowed, facilitating subsequent data-matching.In this experiment, the gridded fisheries data were expanded to 329 entries per month by the inverse-distance weighting method, totaling 32,571 entries for January and May to December 2003-2013.The spatial resolution of the fisheries data was improved from 5 • × 5 • to 1 • × 1 • , providing a more detailed portrayal of the monthly habitat distribution of Bigeye Tuna.The increase in the number of samples reduced the likelihood of random errors in model training, while the gap in spatial resolution with product-level environmental data was narrowed, facilitating subsequent data-matching.

Remote-Sensing Data Processing
The two remote-sensing data used in this experiment need to be preprocessed separately due to the significant differences in data structure and spatiotemporal resolution.
The initially acquired data files for the five marine environmental factors include four attribute dimensions, i.e., time, depth, latitude, and longitude.Eleven marine environmental factors at various depths are generated by reading Env data at specific fixed depth points.The Random Forest algorithm is used in this experiment to fill in a small number of consecutive missing data values in the selected geospatial area due to the presence of missing data in the original Env data.The partial rows and columns without missing values are selected as training data according to the correlation of continuous data.After generating the Random Forest model, the rows and columns with the smallest number of missing values are supplemented first.The dataset is updated after each filling, and this process is repeated until all missing values are filled.
(2) The preprocessing of L1B data The original multispectral remote-sensing data used in this experiment are MODIS L1B level data containing 36 bands (L1B data).Bands 8~16 are related to the remote-sensing information of ocean water color and phytoplankton, which are often used to invert the ocean products of MODIS data, and bands 31~32 are related to the Earth's surface, cloudtop temperatures, and the atmosphere, which can identify invalid values generated by the occlusion of cloudy features.Therefore, the above 11 bands are selected as experimental data to study the relationship between marine remote-sensing data and the distribution of habitats.
The L1B data need to be processed to identify the cloud feature occlusion regions before use to remove invalid factors of non-marine information, and the method selected for this experiment is the Normalized Detection Index method [22].There is a difference between the reflectance of the ocean surface and cloud features in the visible (0.66 µm) and thermal infrared (11 µm) bands.The probability that a point in the ocean is greater than the difference between the two is greater after normalization since the ocean surface has a lower reflectance at 0.66 µm and a higher reflectance at 11 µm.The probability of a point being an ocean is calculated using the Normalized Detection Index method: where P denotes the probability that a point is an ocean, ρ 0.66µm and ρ 11µm denote the reflectivity of the point in the 0.66 µm band and the 11 µm band, and f(x) is a normalization processing function denoting the normalization of x to the interval [-1, 1].A point is considered to be an ocean when it has P ≥ 80%.Several data still exist at the same pixel point each month after identifying the filtered invalid values because the temporal resolution of the L1B data is measured in days.The remaining L1B data for each month were divided into groups of equal days in chronological order in order to harmonize with the temporal resolution of the fishery data, and a mean value was obtained for each group, ultimately retaining 5 data per month to match the fishery data The error introduced by the missing value treatment is greatly reduced by such a treatment, while the features of the L1B data in continuous time are preserved.
(3) Dataset generation The temporal resolution of both Env and fishery data is in months without any further processing.The Env dataset has a spatial resolution of 0.25 • × 0.25 • , while the fisheries data have a spatial resolution of 1 • × 1 • after grid expansion.To make the two datasets compatible, the matrix matching method was used: the Env data were organized into 10×10 size matrices centered on the fisheries' data points [23,24], encompassing a 2.5 • × 2.5 • range of marine environmental conditions around each fisheries' data center point.The dataset generated after matching is referred to as Env_Fishery, with its data structure illustrated in Figure 3.
J. Mar.Sci.Eng.2024, 12, x FOR PEER REVIEW (3) Dataset generation The temporal resolution of both Env and fishery data is in months without ther processing.The Env dataset has a spatial resolution of 0.25° × 0.25°, while the data have a spatial resolution of 1° × 1° after grid expansion.To make the two compatible, the matrix matching method was used: the Env data were organi 10×10 size matrices centered on the fisheries' data points [23,24], encompassing 2.5° range of marine environmental conditions around each fisheries' data cente The dataset generated after matching is referred to as Env_Fishery, with its data s illustrated in Figure 3.The method preserves the original resolution of Env data, effectively reducin lation errors relative to the traditional mean processing matching method.The characteristics of various environmental factors can be accurately reflected Env_Fishery dataset, which is the smallest unit after matching the environment a ery data, providing rich and comprehensive information on environment-fishe ages.
The Env_Fishery dataset has been further organized to temporal patterns in tribution of habitats and the correlation between temporal characteristics of mari ronmental factors and the formation of fisheries.
Studies have shown that the effects of marine environmental factors on the d tion of habitats exhibit a lag, with the impacts of water temperature and phytop content on CPUE values in the current month persisting for several months [25].H the distribution of fishing habitats at the monthly granularity demonstrates tempo tinuity, showing similarities between the same month in different years.Therefor object of study is the distribution of fishing grounds in month m of year y: First, a time series of multiple months within the same year is constructed wi size of 3 months.The Env_Fishery data of m-January, m-February, and m-March y are aligned to correspond with the CPUE value of month m of year y, allowing exploration of the time-series characteristics of consecutive months within a fix Second, a time series of consecutive years with the same month is constructed, se year step to 3 years.The Env_Fishery data of month m in years y−1, y−2, and arranged to correspond to the CPUE values of month m in year y, enabling the exp of the temporal characteristics of consecutive years with fixed months.Finally, poral-environmental-fisheries dataset Env_Fishery_Time is generated, and its da ture is shown in Figure 4.
The Env_Fishery_Time dataset continues to be organized with the Env_Fish The method preserves the original resolution of Env data, effectively reducing calculation errors relative to the traditional mean processing matching method.The spatial characteristics of various environmental factors can be accurately reflected by the Env_Fishery dataset, which is the smallest unit after matching the environment and fishery data, providing rich and comprehensive information on environment-fishery linkages.
The Env_Fishery dataset has been further organized to temporal patterns in the distribution of habitats and the correlation between temporal characteristics of marine environmental factors and the formation of fisheries.
Studies have shown that the effects of marine environmental factors on the distribution of habitats exhibit a lag, with the impacts of water temperature and phytoplankton content on CPUE values in the current month persisting for several months [25].However, the distribution of fishing habitats at the monthly granularity demonstrates temporal continuity, showing similarities between the same month in different years.Therefore, if the object of study is the distribution of fishing grounds in month m of year y: First, a time series of multiple months within the same year is constructed with a step size of 3 months.The Env_Fishery data of m-January, m-February, and m-March of year y are aligned to correspond with the CPUE value of month m of year y, allowing for the exploration of the time-series characteristics of consecutive months within a fixed year.Second, a time series of consecutive years with the same month is constructed, setting the year step to 3 years.The Env_Fishery data of month m in years y−1, y−2, and y−3 are arranged to correspond to the CPUE values of month m in year y, enabling the exploration of the temporal characteristics of consecutive years with fixed months.Finally, the temporalenvironmental-fisheries dataset Env_Fishery_Time is generated, and its data structure is shown in Figure 4.
August to December from 2006 to 2013 as the primary periods for studying fishery habitats due to the requirements of time-series steps in the dataset, and the remaining months serve as supplementary material to explore the relationship between various months and fishery habitats.The five data lines retained each month after L1B data preprocessing are time-continuous and can reflect the dynamics of the MODIS data for that month.Therefore, the step size of the time series was set to 5 to match the fisheries' data information for the current month and harmonize the temporal resolution of the L1B data and the fisheries' data.Spatially, due to the large volume of MODIS data, the use of the matrix matching method may lead to model computational overload.To address this, the averaging method was chosen to adapt to the spatial resolution of the fisheries data and solve the problem of missing invalid values for some L1B data points.Following the matching process, the time-series-multi-spectrum dataset L1B_Fishery_Time is generated, spanning January 2006-2013 and August-December.Its data structure is shown in Figure 5.

Standardization and Normalization of the Experiment Data
The above types of feature data are standardized to have a mean of 0 and a variance of 1 in this experiment in order to unify the magnitude and order of magnitude of the spatiotemporal factors, the environmental factors of the Env data, and the values of each band of the L1B data.The standardized formula is: The Env_Fishery_Time dataset continues to be organized with the Env_Fishery dataset as the base unit and contains feature information in the three dimensions of environment, space, and time.The generated Env_Fishery_Time dataset focuses on January and August to December from 2006 to 2013 as the primary periods for studying fishery habitats due to the requirements of time-series steps in the dataset, and the remaining months serve as supplementary material to explore the relationship between various months and fishery habitats.
The five data lines retained each month after L1B data preprocessing are time-continuous and can reflect the dynamics of the MODIS data for that month.Therefore, the step size of the time series was set to 5 to match the fisheries' data information for the current month and harmonize the temporal resolution of the L1B data and the fisheries' data.Spatially, due to the large volume of MODIS data, the use of the matrix matching method may lead to model computational overload.To address this, the averaging method was chosen to adapt to the spatial resolution of the fisheries data and solve the problem of missing invalid values for some L1B data points.Following the matching process, the time-seriesmulti-spectrum dataset L1B_Fishery_Time is generated, spanning January 2006-2013 and August-December.Its data structure is shown in Figure 5.
August to December from 2006 to 2013 as the primary periods for studying fishery habitats due to the requirements of time-series steps in the dataset, and the remaining months serve as supplementary material to explore the relationship between various months and fishery habitats.The five data lines retained each month after L1B data preprocessing are time-continuous and can reflect the dynamics of the MODIS data for that month.Therefore, the step size of the time series was set to 5 to match the fisheries' data information for the current month and harmonize the temporal resolution of the L1B data and the fisheries' data.Spatially, due to the large volume of MODIS data, the use of the matrix matching method may lead to model computational overload.To address this, the averaging method was chosen to adapt to the spatial resolution of the fisheries data and solve the problem of missing invalid values for some L1B data points.Following the matching process, the time-series-multi-spectrum dataset L1B_Fishery_Time is generated, spanning January 2006-2013 and August-December.Its data structure is shown in Figure 5.

Standardization and Normalization of the Experiment Data
The above types of feature data are standardized to have a mean of 0 and a variance of 1 in this experiment in order to unify the magnitude and order of magnitude of the spatiotemporal factors, the environmental factors of the Env data, and the values of each band of the L1B data.The standardized formula is:

Standardization and Normalization of the Experiment Data
The above types of feature data are standardized to have a mean of 0 and a variance of 1 in this experiment in order to unify the magnitude and order of magnitude of the spatiotemporal factors, the environmental factors of the Env data, and the values of each band of the L1B data.The standardized formula is: where x * denotes the result after standardization of the feature data, x denotes the unstandardized feature data, and µ and σ denote the mean and standard deviation of the type of feature data, respectively.The CPUE values were normalized to the interval [0, 1] in this experiment due to the uneven order-of-magnitude distribution of the CPUE values in the fishery data in order to reduce the error caused by the order-of-magnitude gap.The normalized formula is: y − y min y max − y min (6) where y * denotes the value after CPUE normalization; y denotes the value before CPUE normalization; and y min and y max denote the minimum and maximum values of CPUE in the fishery dataset, respectively.

Habitat Prediction Methods
The temporal characteristics and spatial distribution characteristics of the marine environmental and fishery data in both annual and monthly time dimensions are contained in the Env_Fishery_Time dataset.In addition, the remote-sensing information of the ocean spectrum and its dynamically changing characteristics are provided in the L1B_Fishery_Time dataset as the more original multispectral remote-sensing data.Due to differences in data structure and organization, feature extraction is performed separately on the Env_Fishery_Time and the L1B_Fishery_Time datasets.Decision-level feature fusion is then conducted.The different modules of this experimental model are described below.

Spatiotemporal Feature Extraction Model for Multi-Branch Marine Environmental Information
The spatiotemporal factor, marine environmental factor (Env data), and fishery information (CPUE) are included in the Env_Fishery_Time dataset.The structure of these three parts differs: the spatiotemporal factor data are one-dimensional, while the Env data and fisheries data are multidimensional continuous-time temporal data (consecutive months in the same year and consecutive years in the same month).Therefore, these data need to be input into the model separately and for feature extraction.The structure of the multidimensional spatiotemporal ocean information feature extraction fusion model with the Env_Fishery_Time dataset as input is shown in Figure 6, which is hereinafter referred to as the BiLSTM-1D-Env (Bilayer LSTM-1DCNN-Env) model.
where  * denotes the result after standardization of the feature data,  denotes the unstandardized feature data, and  and σ denote the mean and standard deviation of the type of feature data, respectively.
The CPUE values were normalized to the interval [0, 1] in this experiment due to the uneven order-of-magnitude distribution of the CPUE values in the fishery data in order to reduce the error caused by the order-of-magnitude gap.The normalized formula is: where  * denotes the value after CPUE normalization;  denotes the value before CPUE normalization; and  min and  max denote the minimum and maximum values of CPUE in the fishery dataset, respectively.

Habitat Prediction Methods
The temporal characteristics and spatial distribution characteristics of the marine environmental and fishery data in both annual and monthly time dimensions are contained in the Env_Fishery_Time dataset.In addition, the remote-sensing information of the ocean spectrum and its dynamically changing characteristics are provided in the L1B_Fish-ery_Time dataset as the more original multispectral remote-sensing data.Due to differences in data structure and organization, feature extraction is performed separately on the Env_Fishery_Time and the L1B_Fishery_Time datasets.Decision-level feature fusion is then conducted.The different modules of this experimental model are described below.

Spatiotemporal Feature Extraction Model for Multi-Branch Marine Environmental Information
The spatiotemporal factor, marine environmental factor (Env data), and fishery information (CPUE) are included in the Env_Fishery_Time dataset.The structure of these three parts differs: the spatiotemporal factor data are one-dimensional, while the Env data and fisheries data are multidimensional continuous -time temporal data (consecutive months in the same year and consecutive years in the same month).Therefore, these data need to be input into the model separately and for feature extraction.The structure of the multidimensional spatiotemporal ocean information feature extraction fusion model with the Env_Fishery_Time dataset as input is shown in Figure 6, which is hereinafter referred to as the BiLSTM-1D-Env (Bilayer LSTM-1DCNN-Env) model.The overall strategy of the model involves using a Long Short-Term Memory Network (LSTM) to process multidimensional spatiotemporal data and 1D Convolutional Neural Networks (1D-CNN) for feature extraction on one-dimensional spatiotemporal factors.The marine environmental factor data are spatially represented as a two-dimensional matrix, so the time series feature extraction of high-dimensional data is carried out using the ConvLSTM2D layer.This layer is a convolutional LSTM, similar to the LSTM layer, but with the inputs and transforms in the form of convolutions, allowing it to fully mine the marine information features and their spatial correlations in the two-dimensional space of longitude and latitude.Fishery data are time-series data without spatial structure and can be used directly for time-series feature extraction using the LSTM layer.Both the ConvLSTM2D layer and the LSTM layer are double-layer structures, and the problem of gradient vanishing and gradient explosion caused by single-layer structures to avoid the problem of gradient vanishing and gradient explosion caused by single-layer structures, making it easier to capture the long-term dependence of long-time sequences.The spatiotemporal factor part is chosen to be extracted using 1D-CNN because it is 1-D data and 1D-CNN can effectively extract local features and features of different scales in the data.
There are corresponding feature quantities for spatiotemporal factors, marine environmental factors (Env data), and fishery information (CPUE) after feature extraction.The fully connected network was used to perform feature fusion on the spliced three-part features to obtain the prediction results of CPUE values based on the Env_Fishery_Time dataset.
This module was used to forecast the CPUE distribution for the coming months based on the horizontal and vertical multidimensional time series, as well as the marine environmental factors and fishery information in the surrounding spatial locations.

Feature Extraction Model for Two-Branch Optical Remote-Sensing Data
The L1B_Fishery_Time dataset is composed of two parts: multispectral data (L1B data) and spatiotemporal factors.The structure of the feature extraction fusion model for multi-source remote-sensing data with the L1B_Fishery_Time dataset as the input is shown in Figure 7, which is hereafter referred to as BiLSTM-1D-L1B (Bilayer LSTM-1DCNN-L1B) model.The two parts of feature data used in this experiment (Env_Fishery_Time dataset, L1B_Fishery_Fishery_Time dataset) are quite different in structure and dimension, making it more complicated to carry out feature-level fusion directly with a single model.The decision-level fusion method is used in this experiment to fuse the features of the two parts of the data in order to avoid the loss of effective features and also to prevent the generation of overfitting.The structure of the decision-level fusion model for feature extraction of heterogeneous data with the Env_Fishery_Time and L1B_Fishery_Time datasets as inputs is shown in Figure 8, which is hereinafter referred to as the BiLSTM-DF (Decision Fusion Model) model.Feature extraction model for two-branch optical remote-sensing data with L1B_Fishery_Time dataset as input.
The spatiotemporal factor part is consistent with the Env_Fishery_Time dataset and is processed using the same module.For the multispectral data, two layers of LSTM were used to extract the pattern of change in the distribution of the ocean remote-sensing information as well as the time series information since the multispectral data have 11 bands and the time dimension includes five data points per month, and then the effective sequence features were further extracted using the Dense layer of multilayer neural network.
The corresponding feature quantities are obtained for the spatiotemporal factor and the multispectral data (L1B data) after feature extraction, similarly to the process described in Section 2.3.1.The feature fusion of these two parts is performed using a fully connected network to obtain the prediction results of the CPUE values on the L1B_Fishery_Time dataset.
This module captures the real-time marine environment through multispectral remotesensing information and enriches the real-time features affecting the CPUE distribution of the month.
The two parts of feature data used in this experiment (Env_Fishery_Time dataset, L1B_Fishery_Fishery_Time dataset) are quite different in structure and dimension, making it more complicated to carry out feature-level fusion directly with a single model.The decision-level fusion method is used in this experiment to fuse the features of the two parts of the data in order to avoid the loss of effective features and also to prevent the generation of overfitting.The structure of the decision-level fusion model for feature extraction of heterogeneous data with the Env_Fishery_Time and L1B_Fishery_Time datasets as inputs is shown in Figure 8, which is hereinafter referred to as the BiLSTM-DF (Decision Fusion Model) model.The two parts of feature data used in this experiment (Env_Fishery_Time dataset, L1B_Fishery_Fishery_Time dataset) are quite different in structure and dimension, making it more complicated to carry out feature-level fusion directly with a single model.The decision-level fusion method is used in this experiment to fuse the features of the two parts of the data in order to avoid the loss of effective features and also to prevent the generation of overfitting.The structure of the decision-level fusion model for feature extraction of heterogeneous data with the Env_Fishery_Time and L1B_Fishery_Time datasets as inputs is shown in Figure 8, which is hereinafter referred to as the BiLSTM-DF (Decision Fusion Model) model.The decision-level fusion method used in this module is the weighted average method, calculated as follows: where Y ˆ denotes the final prediction result; Y ˆL1b_Fishery _Time and Y ˆEnv_Fishery_Time denote the prediction results of the fusion model of feature extraction for the two datasets, respectively; ω 1 _ℎ _ and ω  _ℎ _ denote the weights of the prediction results in the final prediction results of the two datasets.Results in the final prediction results and the specific weight values were determined through multiple experiments.
Integration and aggregation of decisions are emphasized by decision -level fusion models rather than merely the integration of features, making them more interpretable relative to single-feature fusion models.These models are also more flexible in dealing with data of different structures and characteristics, which aligns well with the characteristics of the data in this experiment.The decision-level fusion model has shown better results relative to feature-level fusion in this experiment.The decision-level fusion method used in this module is the weighted average method, calculated as follows: where Ŷ denotes the final prediction result; ŶL1b_Fishery_Time and ŶEnv_Fishery_Time denote the prediction results of the fusion model of feature extraction for the two datasets, respectively; ω L1b_Fishery_Time and ω Env_Fishery_Time denote the weights of the prediction results in the final prediction results of the two datasets.Results in the final prediction results and the specific weight values were determined through multiple experiments.
Integration and aggregation of decisions are emphasized by decision-level fusion models rather than merely the integration of features, making them more interpretable relative to single-feature fusion models.These models are also more flexible in dealing with data of different structures and characteristics, which aligns well with the characteristics of the data in this experiment.The decision-level fusion model has shown better results relative to feature-level fusion in this experiment.

Experiment Environment
Two datasets were used in the above experimental setting: the Env_Fishery_Time and L1B_Fishery_Time datasets.The Env_Fishery_Time dataset includes data for January and May to December 2006-2012 and January and May to November 2013, totaling 15,463 records.These data were divided into training and validation sets in the ratio of 3:7, while the test set comprised the data from December 2013 with 329 entries.The L1B_Fishery_Time dataset was segmented similarly.The training and validation sets were used to validate the model's performance, the training and test sets were used to evaluate the generalization performance of the model, and the combined results were used to evaluate and analyze the experiment in general.
The following hyperparameters were selected as the best performing after several comparison attempts: the maximum number of iterations of the model is set to 500 by observing the loss value descent curve to make the loss value of the model converge stably; the learning rate is adjusted using the Adam optimization algorithm so that it can be adjusted adaptively at different stages of training to better adapt to the characteristics of the data; and the value of the Dropout function is set to 0.5 to prevent overfitting.

Model Performance Evaluation Index
The Root Mean Square Error (RMSE) and Goodness of Fit (R 2 ) are used in this experiment to evaluate the model performance.RMSE is used to calculate the mean of the model's prediction error to measure the difference between the predicted and true values.In general, the smaller value of RMSE indicates a smaller difference between the model's prediction and the actual value and thus, better prediction accuracy.;R 2 is used to indicate the goodness of fit and takes values in the range of [0, 1].The closer the R 2 is to 1, the better the model fit.RMSE and R 2 are calculated as: where n denotes the number of samples, ȳ denotes the sample mean, y i and ȳ denote the true and predicted values of the i-th sample, respectively.

The Module of Feature Extraction
The spatial scale of the environmental factor data matrix in the Env_Fishery_Time dataset is analyzed in this experiment, and the LSTM sequence length parameter of the feature extraction model using the Env_Fishery_Time dataset as input is adjusted.
The environmental factor data need to be organized into a data matrix in 2D space to match the CPUE because the spatial resolution of the CPUE in the dataset is in units of 1 • × 1 • , and the environmental factor data are in units of 0.25 • × 0.25 • .So, the environmental factors were organized according to 4 × 4, 6 × 6, 8 × 8, 10 × 10, and 12 × 12, i.e., environmental factor data matrices at different spatial scales of 1  1.
It can be seen from Table 1 that a 1 • × 1 • unit of fishing area can just be covered by a matrix of environmental factor data organized by 4 × 4 in units of 0.25 • × 0.25 • .As the scale of the matrix increases, therefore encompassing a broader range of surrounding marine environments, the prediction error RMSE of the model gradually decreases while its fitting ability improves.This trend indicates that the distribution and changes in the surrounding marine environmental factors significantly influence the CPUE values at the central study site.  1 reaches an optimum of 0.95 or higher when the data matrix is organized as 10 × 10, and RMSE achieves a minimum value of 0.03210, indicating that correlation is strongest for this spatial organization.However, when the data matrix is expanded to 12 × 12, there is a slight decrease in the model's predictive effectiveness, indicating that further increasing the spatial scale does not improve the model's performance.The correlation between the marine environment of the surrounding range and the CPUE values of the central study site weakened as the distance increased.Additionally, the ability of the model's feature extraction and fitting prediction will be affected by the overly complex data structure that tends to overfit the model during feature extraction.
Therefore, the 0.25 • units of 10 × 10 organization of the environmental factor data matrix in the Env_Fishery_Time dataset was used in this experiment to extract the characteristics of the marine environment in and around the study site.
Next, the LSTM sequence length parameter of the feature extraction model with the Env_Fishery_Time dataset as input is tuned.The sequence length is a crucial hyperparameter that impacts the performance and training efficiency of the model when using LSTM for time series prediction.The length of the year sequence is set to 1, 2, and 3 years, and the length of the month sequence is set to 1 month, 2 months, and 3 months.The years are denoted as the first, second, and third years, and the months are denoted as the first, second, and third months for sequences of length 3 from the month of the study.The structure of the multidimensional time series is shown in Figure 9.To maintain balanced feature extraction across different time dimensions, the year and month are extracted with the same sequence length.The years and months are divided into two groups in chronological order: the first year and the first month form the first group, the second year and the second month form the second group, and the third year and the third month form the third group.The correlation between the environmental factors and the CPUE value of the month in each group is analyzed separately, and the results of the correlation analysis are shown in Table 2 (The data in the table are absolute values and mainly used to explore the magnitude of the correlation between CPUE and environmental factors under different subgroups.).The effect of environmental factors on CPUE values has a lag by comparing the data in the table.The correlations between most of the environmental factors (e.g., sst100, Chlα, sla, sss0, sss400, o300) and the CPUE of the current month are strongest in the second group, and the correlations tend to increase and then decrease with increasing time distance.The correlations of a small number of environmental factors (e.g., sst200, sss200, o200) with the CPUE of the month are not significantly different between the first and second groups, but the correlations are significantly lower in the third group.The correlation between a few environmental factors (e.g., sst0) and the CPUE of the current month is larger in the third group than in the first group, but the correlation of the second group is still the strongest, with the overall correlation showing a tendency to increasing and then decreasing.3. It is evident from the comparison that the correlation between the CPUE values of different years and months and the CPUE of the current month diminishes with the increase of the time distance.In other words, the closer the time distance is, the stronger the correlation between the CPUEs is.To further determine the optimal sequence lengths of years and months, comparative experiments were conducted.The results of these experiments, which evaluated the prediction effectiveness of LSTM with different input sequence lengths, are shown in Table 4.It can be seen from Table 4 that the optimal prediction performance of the LSTM model is achieved with a sequence length is 3.The prediction results of the LSTM module are poor due to the fact that there is a lag in the influence of environmental factors on CPUE, and the closer the time distance, the stronger the correlation between CPUE, the effective factors with a strong correlation with the CPUE value of the current month cannot be completely extracted with the length of the time series of one.When the sequence length is increased to 2, the prediction performance improves significantly, likely because the correlation of the environmental factors peaks in the second year and month, which also maintains a certain correlation of the CPUE values.Increasing the sequence length to 3 further enhances prediction results, but the enhancement is significantly reduced, as the correlation between environmental factors and CPUE weakened over longer time distances.

The results of analyzing the correlation between the CPUE values of different years and months and the CPUE values of the current month are shown in Table
Longer sequence lengths allow for better capture of long-term dependencies and trends in the data, which enhances prediction accuracy.However, as the length of the series increases, the organization and structure of the data become more complex, which can create computational cost issues for the model and can also easily lead to overfitting.When the sequence length is set to 4, the RMSE of the model instead increases and R 2 is no longer further improved.The stability of the model also fluctuates due to the increase in computational effort.Therefore, the time sequence length of the LSTM prediction module is set to 3, which performs well both in the training and generalization phases is optimal.

Decision-Level Fusion Module
The weights of the decision-level fusion methods, i.e., the weighted average method, need to be determined experimentally in the decision-level fusion module.There are multiple model-validation experiments and generalization experiments for the BiLSTM-1D-L1B feature extraction model with the L1B_Fishery_Time dataset as input and the BiLSTM-1D-Env feature extraction model with the Env_Fishery_Time dataset as input, respectively, in order to eliminate the chance errors.The results of three randomly selected experiments with the same parameters were recorded, and the results are shown in Table 5.The average values of the different experimental effects are calculated in Table 5 to compare the predictive effectiveness of the two models.In the model-validation experiments, the average prediction error of the BiLSTM-1D-Env model is 0.03317, and the goodness of fit reaches 0.9536.This represents a 12.22% decrease in prediction error compared to 0.03779 of the BiLSTM-1D-L1B model, with the goodness of fit slightly higher than that of the BiLSTM-1D-L1B model, indicating better performance in model validation.In the generalization experiments, the BiLSTM-1D-L1B model exhibits an average error of 0.08807, which is a significant decrease compared to the 0.1016 of the BiLSTM-1D-Env model.The R 2 value for the BiLSTM-1D-L1B model shows a 28.86% improvement, reaching 0.3846 compared to the BiLSTM-1D-Env model.These differences in experimental outcomes between the two heterogeneous data feature extraction models are influenced not only by the models themselves but also by their dependency on different datasets.
The optimal weights for the model-validation experiments and generalization experiments are determined by analyzing the differences in the performance between the feature extraction models for the two datasets.Specifically, the method involves matching the three experimental results of the BiLSTM-1D-Env and BiLSTM-1D-L1B models from Table 5, creating 9 combination methods.The weighting calculation is carried out nine times according to the set weights, and the average performance of the nine feature-level fusion results under the weights is recorded.The average results of decision-level fusion for similar experiments with different weight ratios are shown in Table 6.The weight given to the experimental results of a particular feature extraction model represents the proportion of that model's results in the total results, as well as the proportion of that dataset's features in the total set of feature results.The experiment results weight of the Env_Fishery_Time dataset is assumed to be increasing, and the experimental effect of decision-level fusion is then gradually transitioned from the L1B_Fishery_Time dataset to the Env_Fishery_Time dataset.The experiment result of decision-level fusion represents the average effect of the Env_Fishery_Time dataset and vice versa.Analyzing the experimental results in the table reveals that in this transition process, the experimental error RMSE does not decrease linearly from high to low but decreases initially and then increases.This indicates the presence of an optimal proportional weight that minimizes the error results of the decision-level fusion experiments, making them smaller than the results of the experiments on feature extraction for the two datasets.For instance, the transition process of the generalization experiment for the medium RMSE change curve is shown in Figure 10.The results in Table 5 are analyzed in conjunction with the fact that the Env_Fishery_ Time dataset slightly outperforms the L1B_Fishery_Time dataset in the model-validation experiments-though the difference is not substantial-suggests that decision fusion is optimized when the two weights are approximately equal.It can be seen from Table 6 that the RMSE of the model-validation experiment is minimized at 0.03136 when the weights are balanced, which corroborates this finding.However, the L1B_Fishery_Time dataset performs better in the generalization experiments, showing significantly better results in both the error and the goodness-of-fit metrics compared to the Env_Fishery_Time dataset.The weight of the feature extraction results from the L1B_Fishery_Time dataset needs to be increased to achieve optimal results.It can be seen from Table 6 that the best generalization results are obtained when the weight ratio is set to 6:4, with the RMSE decreasing to 0.07824 and the goodness-of-fit R 2 reaching 0.6374.This represents a 15.17% improvement over the results obtained using the L1B_Fishery_Time dataset alone, thus validating the analysis.The decision-level features fusion effectively integrates features from different datasets, resulting in better experiment outcomes than those obtained from single-dataset feature extraction.Therefore, the optimal weight ratio values for the decision-level fusion module are selected to complete the model-validation experiments and generalization experiments, therefore confirming the effectiveness of the decision-level fusion approach.

Comparative Analysis of the Results of Data-Organization Validation and Model Validation
To validate the data organization and the effectiveness of the models used in this experiment, we compared the results of each model in the model-validation experiment [26][27][28][29].The comparison includes the Bayesian Regression Model (BR), the Regression Tree Model (RT), the Random Forest (RF) model, the Back Propagation Neural Network (BPNN), the 1D Convolutional Neural Network Model (1D-CNN), a sequence model combined from two neural network models (CNN-LSTM), a two-layer time series prediction model with separate inputs for spatiotemporal factors and environmental data (LSTM-LSTM-BPNN), and two feature-level fusion models derived from heterogeneous data.The results of these comparisons are presented in Table 7.The Env_Fishery dataset without time-series features is used for the traditional base models to complete the comparison test, as these models do not support sequence prediction.Dimensionality reduction of the Env_Fishery dataset is required for models that do not support high-dimensional feature datasets.It can be seen from the results in Table 7 that the traditional statistical model represented by BR has a large bias in the prediction of the habitat, with an RMSE of around 0.14.These models struggle to accurately capture the correlation between environmental factors and CPUE.
In the machine-learning model, the RT model uses a recursive division of features to achieve regression prediction of sample values.After adjusting the number of leaf nodes and the height of the tree, the fitting effect is further improved compared to the statistical model, with R 2 increasing from 0.17 to 0.33.The RF model based on the integrated learning of multiple decision trees can further improve the prediction effect.The error decreases significantly compared to the RT model, with the RMSE dropping from 0.12319 to 0.08828.Overall, there is a large improvement in the prediction effect of the machine-learning model compared to the traditional mathematical model.However, the generalization ability of both RT and RF models is poor due to the poor stability of the RT model itself, and the changes in input features can easily lead to overfitting.
The overfitting can be effectively avoided in deep-learning-based neural network models with the introduction of the Dropout layer.The issues arising from the gradient propagation process can also be reduced by the Batch Normalization (BN) layer, which enhances the generalization ability of the model.As shown in Table 7, the R 2 of both neural network models, BPNN and 1D-CNN, can be kept above 0.65, with the 1D-CNN model slightly outperforms BPNN.The reason for this is that 1D-CNN is able to efficiently capture the local patterns in the input sequence data and capture the long-distance dependencies by stacking convolutional layers, and thus 1D-CNN is more suitable for 1D-CNN compared to BPNN for processing sequence data.
Since both CPUE and marine environmental factors are time-series data, in order to further explore the characteristics of their time series, this experiment uses two combined sequential models for comparison: a sequence model combining CNN and LSTM (CNN-LSTM) and a two-layer time-series prediction model with separate inputs for spatiotemporal factors and environmental data (LSTM-LSTM-BPNN).Predictions are made using the Env_Fishery_Time dataset, which contains spatiotemporal features, and the L1B_Fishery_time dataset.The effectiveness of time-series prediction in habitat prediction was demonstrated by the experiments, with the RMSE of both models dropping below 0.058 and R 2 improving to a maximum of 0.9157.The convolutional structure in the CNN-LSTM model can adequately extract Env data features in 2D space, while the LSTM structure in LSTM-LSTM-BPNN is more advantageous for extracting L1B data features that do not contain spatial dimensions.Overall, the LSTM-LSTM-BPNN model, with separate inputs of spatiotemporal and environmental data, shows better performance.Therefore, it is necessary to select appropriate model structures for heterogeneous data to fully mine the data features of each part of the heterogeneous data because a single model cannot meet the requirements of feature extraction for different data structures.
Based on the comparative analysis of the above model experiments, this experiment further improves the LSTM-LSTM-BPNN model and proposes the BiLSTM-DF model to address the problem of heterogeneous data structures with different characteristics.The model includes separate appropriate modules for different parts of the data structure on a single dataset.For the L1B dataset, which benefits from feature extraction via the twolayer LSTM structure, this configuration is retained.For Env data with 2D spatial features, the LSTM layer is replaced with a ConvLSTM2D layer.The ConvLSTM2D captures the temporal and spatial correlations of the Env data by performing convolution operations at each time step to extract the local features in the input sequence.Additionally, the 1D-CNN model, which performs better in processing 1D data processing, is used for feature extraction of spatiotemporal factors in the dataset.Finally, the feature extraction results on both datasets are fused at the decision level using a weighted average algorithm.
In Table 6, the performance of the BiLSTM-DF model on the two datasets (BiLSTM-Env model and BiLSTM-L1B model) is significantly improved compared to the simple sequence model, with the prediction errors decreasing to 0.03317 and 0.03779, respectively.Further comparison of the results of the decision-level fusion (BiLSTM-DF model) of the two datasets reveals that the decision-level fusion method proposed in this experiment is significantly better than the simple feature-level fusion method.The result of feature fusion achieves an R 2 of 0.96278, and the RMSE is reduced to 0.031361, which effectively solves the problem of the difficulty in fusing the features of the data with different dimensions and magnitudes.

Analysis of the Results of the Generalization Experiments
The LSTM-LSTM-BPNN model, which performed well (R 2 > 0.9) on both modelvalidation experiments on both datasets, is selected for comparison experiments in this experiment to validate the generalization ability of the BiLSTM-DF model.The results of the comparison experiments are shown in Table 8.
The data from December 2013 is selected as the test sample of generalization ability in this experiment.Table 8 shows that the proposed model demonstrates superior generalization capability compared to the simple combinatorial sequence model with separate inputs for spatiotemporal factors.Specifically, the RMSE of the BiLSTM-1D-Env model is reduced by 7.1%, and R 2 is increased by 0.083 compared to the LSTM-LSTM-BPNN model.Similarly, the RMSE of the BiLSTM-1D-L1b model decreased by 6.6%, and R 2 increased by 0.065 relative to the LSTM-LSTM-BPNN model.Additionally, the generalization results of the BiLSTM-DF model are also significantly better than the BiLSTM-1D-Env-L1b model by comparing the fusion performance of the two feature extraction model results.The decision-level fusion approach achieves an R 2 of 0.6374 and reduces the RMSE to 0.07824, which indicates that the BiLSTM-DF model proposed in this experiment has a good generalization ability.The model can accurately predict habitat conditions for the upcoming months, with prediction errors not exceeding 0.078.There is a big improvement in this part by the BiLSTM-1D-Env model.Although predicting extremely high values remains challenging, the prediction error for habitats with low CPUE values is greatly reduced, and the prediction curves align more closely with the target curves.In comparison to the scatter plot distributions, the predicted scatter distribution of the BiLSTM-1D-Env model is more centrally converged, especially for CPUE values between 0.5 and 0.8, where the error is reduced compared to the LSTM-LSTM-BPNN model.The prediction of the smaller values is also significantly improved, consistent with the analysis from the prediction curve graph.
The visualization results of the LSTM-LSTM-BPNN and BiLSTM-1D-L1B models on the L1B_Fishery_time dataset are shown in Figure 12 to Figure 13.The prediction curve graphs of both models are analyzed, and the prediction results of CPUE in the LSTM-LSTM-BPNN model are basically in line with the trend of the real value of CPUE, but the overall error is larger, and some samples exhibit significant errors.The improved BiLSTM-1D-L1B model is based on the LSTM-LSTM-BPNN model, and the overall error is reduced.The CPUE prediction curve of the BiLSTM-1D-L1B model aligns more closely with the true value curve, and the number of samples with large errors is significantly reduced.The above phenomenon is also reflected in the scatter plot.The LSTM-LSTM-BPNN model's scatter distribution shows a clustering bias for some samples: the predicted values for samples in the 0.5 to 0.6 range are high, while those for samples in the 0.6 to 0. However, it can still be seen in Figure 15 that the BiLSTM-DF (Decision Fusion Model) model is not effective in predicting the larger value samples, and there are still individual samples with a large gap in the true value.The next step is to further analyze the error of the prediction results from the visualized heat map comparing the prediction results in December 2013 with the real distribution of the habitat.
The comparison results of the visualized heat map for December 2013, shown in Figures 16 and 17, reveal that the predicted habitat distribution generally aligns with the real fishery distribution, clearly displaying both central and non-central habitats.In the real distribution map of the fishing habitats, it was found that there are two main center fishing habitats at 5 • W, 45 • S and 5 • W, 80 • S, which are noted as Habitat 1 and Habitat 2, respectively; There is one sub-center of the habitat, at 5 • W, 65 • S, noted as Habitat 3. It can be found by analyzing the visualized heat map of the prediction results that the locations of the fishing habitats 1 and 3 are basically predicted correctly, and the deviation of the location is within 1 • .The prediction of Habitat 2 was more ambiguous, and although it was possible to predict the location of a wide range of habitats, it was not possible to make an accurate prediction of the center point of the Habitat.It can be observed by analyzing the values of the prediction results that for the larger values in the central habitat, the prediction results are small, which is consistent with the results obtained from the analysis in Figure 15.
cation is within 1°.The prediction of Habitat 2 was more ambiguous, and although it was possible to predict the location of a wide range of habitats, it was not possible to make an accurate prediction of the center point of the Habitat.It can be observed by analyzing the values of the prediction results that for the larger values in the central habitat, the prediction results are small, which is consistent with the results obtained from the analysis in Figure 15.The following are possible reasons for the error based on the prediction performance in the generalization experiments: (1) There are some limitations in the correlation between the two remote-sensing datasets and the CPUE of the habitat.While marine environment factors and marine remote-sensing information are mainly used to predict the distribution of CPUE, better prediction results are obtained due to the many factors affecting the distribution of the habitat.Specifically, the impact of the ocean monsoons and the migratory characteristics of the survival cycle of the fish itself on the accuracy of the habitat prediction is still not ruled out.(2) There might be a chance as the generalization experiment of this experiment only used one month, December 2013, as the test data.

Discussion
As a highly migratory fish, Bigeye Tuna is usually able to sensitively perceive environmental changes in order to determine the direction of migration.The migration of migratory fish is seasonal, and changes in environmental factors under the same or adjacent seasons will have an impact on it.Current research on habitat prediction is mainly realized by exploring the correlation between environmental factors and fishery formation in the current month or using mathematical models to make temporal prediction s of small-scale fisheries based only on changes in the fishery abundance indicator CPUE.As for Bigeye Tuna, as a fish with migratory characteristics and seasonal activity patterns in pelagic fishery, its distribution and migratory behavior are affected by the continuity of many marine The following are possible reasons for the error based on the prediction performance in the generalization experiments: (1) There are some limitations in the correlation between the two remote-sensing datasets and the CPUE of the habitat.While marine environment factors and marine remotesensing information are mainly used to predict the distribution of CPUE, better prediction results are obtained due to the many factors affecting the distribution of the habitat.Specifically, the impact of the ocean monsoons and the migratory characteristics of the survival cycle of the fish itself on the accuracy of the habitat prediction is still not ruled out.(2) There might be a chance as the generalization experiment of this experiment only used one month, December 2013, as the test data.

Discussion
As a highly migratory fish, Bigeye Tuna is usually able to sensitively perceive environmental changes in order to determine the direction of migration.The migration of migratory fish is seasonal, and changes in environmental factors under the same or adjacent seasons will have an impact on it.Current research on habitat prediction is mainly realized by exploring the correlation between environmental factors and fishery formation in the current month or using mathematical models to make temporal predictions of small-scale fisheries based only on changes in the fishery abundance indicator CPUE.As for Bigeye Tuna, as a fish with migratory characteristics and seasonal activity patterns in pelagic fishery, its distribution and migratory behavior are affected by the continuity of many marine environmental factors, and the small-scale habitat prediction of the fishery does not meet the needs of the pelagic fishery.Meanwhile, in terms of the use of marine environmental data and feature extraction methods, although the sequence model feature extraction methods with single-source data as input can reflect the changes in the marine environment to a certain extent, it is often difficult to capture the deep mechanisms behind the complex and variable fish migration behavior.
To explore the correlation between ocean remote-sensing information and Bigeye Tuna habitat distribution in a multidimensional spatiotemporal environment and to better predict its habitat prediction, this paper constructs time-series heterogeneous remotesensing datasets based on the environmental remote-sensing data and the original optical remote-sensing data, respectively, and proposes a time-series prediction model that can extract, fuse, and predict spatiotemporal features of the heterogeneous data in a multidimensional way, i.e., the BiLSTM-DF (Decision Fusion Model) model.The model is based on the LSTM time-series structure for deep feature extraction in four dimensions of time-space-environment-spectrum (TSES) for two heterogeneous sources of remotesensing data, and the dynamic spatiotemporal features affecting the distribution of fishing grounds are mined in depth.On the product-level dataset of marine environmental factors, a multi-branch spatiotemporal feature extraction model of marine environmental information based on product-level remote-sensing data are constructed, and the time-series features and two-dimensional spatial features of the marine environmental factor part of the dataset are extracted using the ConvLSTM2D structure; Mining of temporal variations of CPUE in previous years and months using LSTM layer; extraction of spatiotemporal factor features of fishing grounds using 1D-CNN structure; In the multi-source remotesensing dataset, a two-branch optical remote-sensing data feature extraction model is constructed to extract the spectral time-series features of the daily changes in the marine environment, and the LSTM model is used to explain the effects of the dynamic changes in the information of different marine bands on the CPUE on the time series; Finally, the results of the two data feature extraction are fused at the decision level, and the optimal weight values are determined experimentally using a weighted average method to obtain CPUE predictions for future months.In the model-validation experiments, the R2 of the BiLSTM-DF (Decision Fusion Model) model reaches 0.96278, and the prediction error is reduced to 0.031361, which has a significant advantage over the simple feature fusion model and the combined temporal model; in the generalization experiments, its prediction error is reduced to 0.07824, and the R2 is improved to 0.6374, which can provide accurate habitat prediction for future months.

Conclusions
To address the lack of spatiotemporal dimension correlation studies, the lack of studies on future habitat prediction, the single-ocean remote-sensing dataset, and the difficulty of heterogeneous data fusion in the current fishery prediction, and also to explore the correlation between ocean remote-sensing information and Bigeye Tuna habitat distribution in a multidimensional spatiotemporal environment, and to realize a more accurate habitat prediction for Bigeye Tuna, we take the Bigeye Tuna fishery data of the Southwest Indian Ocean (SWIO) as the material for fishery prediction studies based on the multi-feature fusion of heterogenous remote-sensing data.The Bigeye Tuna fishery data in the Southwest Indian Ocean is used as the material for the research of habitat prediction based on multifeature fusion of remote-sensing data from heterogeneous sources.This paper breaks through the limitation of single-source data and proposes to utilize heterosource remotesensing data for time-space-environment-spectrum (TSES) four-dimensional deep feature extraction, and realizes the deep mining of dynamic spatiotemporal features affecting the distribution of fishery by constructing a multi-branch temporal feature extraction model, and the specific conclusions of the study are as follows: 1.
Aiming at the lack of research on the correlation of spatiotemporal dimensions and the lack of research on the prediction of future habitat prediction, this paper proposes a method of habitat prediction based on multi-branch spatiotemporal feature extraction of marine environmental information from environmental remote-sensing data, which provides a deep excavation of dynamic spatiotemporal features influencing the distribution of the habitat and is able to provide a distributional prediction of the habitat situation in the future period, and generally obtains good results, and provides a basis for predicting the fishing habitat of tuna species and other migratory characteristic fish species and other migratory characteristics of fish species habitat prediction.

2.
Aiming at the problem of a single marine remote-sensing dataset and the difficulty of heterogeneous data fusion, this paper further joins the MODIS L1B level data and proposes a two-branch optical remote-sensing data feature extraction model based on the L1B level data, which continues to supplement the spatiotemporal dynamic features of the marine environment; it is also proposed to use the multi-feature fusion method of heterogeneous remote-sensing data to fuse the feature extraction results of product-level remote-sensing data and L1B-level remote-sensing data at the decision-making level, which improves the prediction effect of future habitat prediction remarkably, and also effectively solves the problem of not being able to effectively fuse multi-source heterogeneous data due to the differences in the data structure and data scale.

3.
Overall, the prediction effect R 2 of the method proposed in this paper reaches 0.96278, and the RMSE is reduced to 0.031361, which is a large degree of improvement compared with both traditional machine-learning methods and deep learning methods, proving the feasibility, scientificity, and advancement of the fishery habitat-temporal prediction based on the fusion of heterogeneous data, and realizing the in-depth excavation of the distribution factors of fisheries under the influence of multidimensional spatiotemporal environment.It provides a new idea for the habitat prediction of tuna species and other migratory characteristic species.
Finally, there are still some areas of improvement that can be made to this paper.First of all, due to the many factors affecting the distribution of the habitat prediction, this experiment mainly uses marine environmental factors and marine remote-sensing information to predict the distribution of the CPUE, and although better prediction results were obtained, it still does not exclude that the marine monsoon and the migratory characteristics of the survival cycle of the fish itself, etc. have an impact on the accuracy of the habitat prediction.The biological characteristics of fish can be further investigated in subsequent studies to improve the prediction accuracy of the model in generalization experiments.Second, since the generalization experiments in this experiment only use one month in December 2013 as test data, there may be a chance that the scope of the generalization experiments can be expanded to further study the pattern and optimization of habitat prediction.

Figure 1 .
Figure 1.The Bigeye Tuna single-month habitat data before grid expansion in January 2006.

Figure 2 .
Figure 2. The Bigeye Tuna single-month habitat data after grid expansion in January 2006.

Figure 1 .
Figure 1.The Bigeye Tuna single-month habitat data before grid expansion in January 2006.

Figure 1 .
Figure 1.The Bigeye Tuna single-month habitat data before grid expansion in January 2006.

Figure 2 .
Figure 2. The Bigeye Tuna single-month habitat data after grid expansion in January 2006.

Figure 2 .
Figure 2. The Bigeye Tuna single-month habitat data after grid expansion in January 2006.

Figure 6 .
Figure 6.A multi-branch spatiotemporal feature extraction model for marine environmental information using the Env_Fishery_Time dataset as an input.

Figure 6 .
Figure 6.A multi-branch spatiotemporal feature extraction model for marine environmental information using the Env_Fishery_Time dataset as an input.

Figure 7 .
Figure 7. Feature extraction model for two-branch optical remote-sensing data with L1B_Fish-ery_Time dataset as input.

Figure 7 .
Figure 7.Feature extraction model for two-branch optical remote-sensing data with L1B_Fishery_Time dataset as input.

Figure 7 .
Figure 7. Feature extraction model for two-branch optical remote-sensing data with L1B_Fish-ery_Time dataset as input.

Figure 8 .
Figure 8. Decision-level fusion model for heterogenous remote-sensing data.

Figure 8 .
Figure 8. Decision-level fusion model for heterogenous remote-sensing data.

3 . 1 .
Experiment Process 3.1.1.Experiment Environment Windows 10 operating system is used in this experimental workstation, and the TensorFlow 2.1.0framework environment based on Python 3.6 has been set up on the Anaconda3 platform to conduct the experiments with the built-in Keras library with version 2.3.1.

5 •
, and 3 • × 3 • were organized to conduct comparison experiments.The comparative experimental results of ConvLSTM2D feature extraction at different spatial scale organization of environmental factors are shown in Table

Figure 9 .
Figure 9. Structural diagram of multidimensional time series.

Figure 9 .
Figure 9. Structural diagram of multidimensional time series.

Figure 10 .
Figure 10.RMSE variation curves under different weight ratios in the generalization experiment.

Figure 15 .
Figure 15.(a) Prediction curves for BiLSTM-DF (Decision Fusion Model) model generalization experiments; (b) Scatter plots for BiLSTM-DF (Decision Fusion Model) model generalization experiments.The visualization results of the LSTM-LSTM-BPNN and BiLSTM-1D-Env models on the Env_Fishery_time dataset are shown in Figures 11 and 12.It can be seen from the prediction curve graph that the LSTM-LSTM-BPNN model predicts CPUE results primarily in the range of 0.5 to 0.8, showing poor prediction ability for values outside this interval.There is a big improvement in this part by the BiLSTM-1D-Env model.Although predicting extremely high values remains challenging, the prediction error for habitats with low CPUE values is greatly reduced, and the prediction curves align more closely with the target curves.In comparison to the scatter plot distributions, the predicted scatter distribution of the BiLSTM-1D-Env model is more centrally converged, especially for CPUE values 7 range are low.The overall error problem is solved by the improved BiLSTM-1D-L1B model, with prediction results in the 0.5 to 0.7 range being consistent with the true values and the larger error sample points showing convergence.The results of the decision-level fusion BiLSTM-DF (Decision Fusion Model) model visualization of the two models are shown in Figure 15.The limitations of each of the two datasets and feature extraction models are addressed by decision-level fusion, and the shortcomings of the respective results are complemented by the proportional fusion of the two results.In the prediction curves, the BiLSTM-DF model shows improved prediction accuracy for samples with smaller values, mitigating the excessive errors observed in both individual feature extraction models.The scatter plot further illustrates this improvement: the predicted scatter points of the BiLSTM-DF model are more closely centered around the function's straight line, indicating a reduction in overall error compared to the BiLSTM-1D-Env model results.The scatter distribution of the BiLSTM-DF model is more clustered compared to the results of the BiLSTM-1D-L1B model, which suggests that the predicted values of the individual samples with larger errors converge further to the true values.The effectiveness of the BiLSTM-DF (Decision Fusion Model) model is demonstrated by both visualization results.

Figure 16 .
Figure 16.Heat map visualizing the distribution of December 2013 habitat prediction.

Figure 17 .
Figure 17.Heat map visualizing the true distribution of fishing habitats in December 2013.

Table 1 .
The results of comparative experiments on different spatial scales of the environmental factors' organization styles.

Table 2 .
The correlation analysis between the environmental factors and the CPUE value of the month in different years and months.

Table 3 .
The correlation analysis between the e CPUE value and the CPUE value of the month in different years and months.

Table 4 .
The comparison of experimental results of different input sequence lengths on the prediction effect of LSTM.

Table 5 .
The results of two heterogeneous data feature extraction models in model-validation experiments and generalization experiments.

Table 6 .
The average results of decision-level fusion for similar experiments with different weight ratios.

Table 7 .
Comparison of experimental results for each model in model-validation experiments.

Table 8 .
The results of each model in the generalization experiment.