dataesr/scientific-paragraphs-categorization huggingface dataset:
Here, we chose k to be the true number of cell populations, which we already know in advance. In real applications, that number may not be known and has to be estimated with prior knowledge or other heuristics. Thus, our heuristic here only works if the practitioner already knows the true number of populations. Additional file 1: Figure  S5  shows the performance all methods on all datasets when parameters are tuned by maximizing the silhouette heuristic, and Table  2  summarizes the mean performance across datasets. The UMAP projection of the resulting DRs can be seen in Fig. 4. We can see that the silhouette heuristic works very well for ZinbWave, where it always identifies parameters equal or very close to the best ones in terms of silhouette, and is comparable to the default parameters in terms of AMI.
It also works well for all methods for simple datasets like sc_10x, sc_10x_5cl and sc_celseq2_5cl, where it also identifies parameters equal or close to the best ones for the silhouette score. As shown by the good AMI performance, these are cases where the initial k-means clustering recovers the correct clustering with good accuracy. However, there are also cases of surprising failures on easy datasets, for example for DCA on TabulaMuris, where the parameter set selected by the silhouette heuristic has a very bad AMI and silhouette with respect to the true labels, probably because the initial clustering selected by the silhouette heuristic completely fails Fig. 4  UMAP representation of Zhengmix8eq after DR by each method (in column) using the ANOVA AMI (top two rows) or empirical silhouette (bottom two rows) heuristic to tune parameters.
In each row, cells are colored either based on their true cell type (rows 1 and 3) or based on a k-means clustering to identify the cell types but nevertheless leads to a good empirical silhouette, while simply using the default parameters gives an almost perfect clustering in terms of AMI and a decent silhouette. On the more challenging immune cell datasets, on the other hand, the silhouette heuristic does not seem to be useful (except for ZinbWave). It leads to worse parameters than the default ones for all methods but ZinbWave on the difficult Zheng-mix8uneq and Zhengmix5eq datasets, except for scran on the later one. For the easier Zhengmix4eq and Zhengmix4uneq, it leads to better parameters for all methods but scVI. In summary, we see that automatically tuning parameters to try to increase the silhouette using the silhouette heuristic only works well on relatively simple cases, up to possible dramatic errors, but on more challenging situations where there is no clear separation between cell types, then it can lead to disastrous choices by overfitting a bad initial clustering.
ZinbWave is an exception where, in our benchmark, the silhouette heuristic gives consistently good results. Proposing other heuristics that really help tune parameters is an important open challenge. In order to stimulate further research in that direction, we provide as supplementary data an original dataset consisting of 5000 embeddings (100 for each single cell dataset and for each method, corresponding to 100 random choices of different parameters, including the default ones and the ones achieving the best AMI and silhouette).
Discussion In this study, we have systematically compared the performances of five representative and popular DR methods over ten datasets with known experimental ground truth, representing various levels of biological complexity. Importantly, we have extensively investigated how the choice of parameters for these methods influence their performances and discussed various ways to properly tune parameters. This can inform practitioners about both the capacity of these methods, as well as on the amount of work required to properly tune them. When properly tuned, we did not observe huge differences in performance between the methods, particularly in terms of AMI. Both PCA-based methods (scran and Seurat) are nevertheless outperformed by ZinbWave and both autoencoder-based methods (DCA and scVI), particularly in terms of silhouette.
An interesting result of our study is the large drop in performance, for all methods, on the immune cell mix data compared to the cell lines. This shows that good performances in the later does not necessarily translate to good performances in the former. In particular, the performances on Zhengmix8eq and Zhengmix8uneq showed that the methods failed to properly separate the various T cells populations, probably due to their relative similarity compared to the other cell populations present in the datasets. Being able to separate similar populations is of utmost relevance when investigating, for example, early tumor development, where the tumor cell population still displays a transcriptome very similar to that of cells of the organ of origin. For example in the case of triple negative breast cancer, in which tumor cells originate from normal luminal cell populations, it is crucial to be able to distinguish the various states the tumor cells undergoes towards full transformation, in order to properly target these cells at an early stage of the disease.
Availability of data and materials The four cell lines datasets come from  [12]  and were downloaded from https://github. com/LuyiTian/sc_mixology in April 2019. The TabulaMuris dataset is an in silico mixture containing all the cells from four tissues sequenced with Smart-Seq2 from the Tabula Muris consortium  [22]  and was downloaded with the TabulaMurisData Bioconductor package, version 1. 2. 0. Zhengmix4eq, Zhengmix4uneq, and Zhengmix8eq come from  [11]  and were downloaded from the DuoClustering2018 Bioconductor package, version 1. 2. 0. We generated Zhengmix5eq and Zhengmix8uneq following the same procedure in order to have more complex datasets. All datasets were subject to the same quality control pipeline, using scater  [16] , version 1. 14. 6. We removed cells three median absolute deviations (MAD) under the mean in counts and expressed genes, as well as those three MAD above the mean in percentage of mitochondrial reads.
noltei-dominated meadows at several sites along the European Atlantic coast. To analyze the seagrass seasonal cycle and to select images corresponding to its maximal development, a two-year S2 dataset was acquired for a French study site in Bourgneuf Bay. The potential of S2 to characterize the Z. noltei seasonal cycle was demonstrated for exposed intertidal meadows. The SPC map that best represented seagrass growth annual maximum was validated using in situ measurements, resulting in a root mean square difference of 14%. The SPC and SB maps displayed a patchy distribution, influenced by emersion time, mudflat topology, and seagrass growth pattern. The ability of S2 to measure the surface area of different classes of seagrass cover was investigated, and surface metrics based on seagrass areas with SPC ≥ 50% and SPC ≥ 80% were computed to estimate the interannual variation in the areal extent of the meadow.
In situ seagrass percent cover and biomass measurements A total of 131 stations ranging from 0 (bare sediment) to 100% seagrass cover from the three study sites were randomly sampled (Table  1 ). At each station, SPC over a 20 cm-diameter circle, corresponding to the core area used for biomass sampling (Water Framework Directive -WFD, 2000/60/EC-) (Fig. 3a ), was measured. For each core surface area, a nadir-viewing photograph was acquired, from which seagrass percent cover (SPC cores ) was computed using the ImageJ software (Diaz-Pulido et al. , 2011). To estimate seagrass leaf biomass (SB cores ), the sediment below the core area was sampled to a depth of approximately 20 cm. Each sample was sieved with a 1 mm mesh and the seagrass leaves were weighted after drying for 48 h at 60 °C  (Bargain et al.
, 2012; Barillé et al. , 2010). Seagrass cover measurements were also performed for satellite retrieval validation along three 1 km transects in Bourgneuf Bay, similar to the sampling designs of previous match-up exercises  (Phinn et al. , 2008; Lyons et al. , 2011; Roelfsema et al. , 2014). For each transect, sampling stations were separated by at least 40 m to avoid autocorrelation between adjacent stations. The stations corresponded to the center of S2 pixels and were geolocated at submetric accuracy using a Trimble® Geo 7 dGPS. For each pixel, five 0. 25 m 2 squares were positioned to form a cross within an area of approximately 25 m 2 so as to take spatial patchiness into account  (Burdick and Kendrick, 2001). For each square, SPC was estimated from a nadir-viewing photograph using the ImageJ software.
In situ radiometry The sky was cloud-free during the acquisition of all radiometric data, and all measurements were performed avoiding shadows over the targets. The upwelling radiance (L core , W m -2 sr -1 nm -1 ) was measured at nadir from 350 to 2500 nm, at the center of the core surface, using an Analyzed Spectral Device field portable spectroradiometer (ASD Fieldspec). The downwelling radiance (L reference , W m -2 sr -1 nm -1 ) diffusively reflected by a Spectralon white reference was also measured. The reflectance of the seagrass cover (R insitu , dimensionless) was estimated following Eq. (1)  (Milton et al. , 2007). = R L L insitu core reference (1) The final dataset consisted of reflectance spectra spanning a range of SPC, from 0 to 100%. As the field work was performed within four-hour periods (corresponding to +/-2 h of low tide) and at different latitudes, the sun elevation and viewing geometry varied within the dataset.
As a result, the raw reflectance spectra exhibited a large range of variability, due either to differences in seagrass cover and/or in acquisition geometry. In order to minimize the uncertainty caused by differences in measurement conditions, a multiplicative scatter correction (MSC;  Isaksson and Kowalski, 1993)  previously applied to seagrass reflectance data  (Fyfe, 2003)  was adapted for application to our in situ datasets (see appendix for more details). Once corrected, the corresponding reflectance spectra were degraded to the spectral resolution of Sentinel-2A using its spectral response function  (SRF; ESA, 2015).
Selection of the vegetation index Many different vegetation indices (VIs) have previously been applied to multispectral remote sensing images to map intertidal seagrass beds  (Bargain et al. , 2012). In the present study, several VIs suitable for S2 were tested using the Bourgneuf Bay datasets. An evaluation was performed for each of the following: normalized difference vegetation index (NDVI;  Tucker, 1979) , normalized difference aquatic vegetation index (NDAVI;  Villa et al. , 2013 Villa et al. , , 2014)) , water adjusted vegetation index (WAVI;  Villa et al. , 2014) , soil-adjusted vegetation index (SAVI;  Huete, 1988) , atmospherically resistant vegetation index (ARVI; Kaufman and Tanre, 1992), modified narrow-band NDVI (mNDVI;  Bargain et al. , 2012) , and modified normalized difference (mND;  Sims and Gamon, 2002).
We evaluated the robustness of each index in terms of differences in the VI vs. SPC cores relationship between 2018 and 2019, using an analysis of covariance (ANCOVA) (Table  A1 ). In all cases, a linear regression was obtained. Only those VIs for which significant differences in the slope and intercept between 2018 and 2019 were not found were selected. We then merged the Bourgneuf Bay 2018 and 2019 datasets, and assessed the performance of the linear regression using the coefficient of determination (R 2 ) and the root mean squared difference (RMSD). The best performance was achieved by the NDVI(665,842), NDVI(705,842), NDAVI(490,842), and ARVI(490,665,842) (Table  A1 ). Besides its good performance, the NDVI(665,842) (hereafter NDVI cores ) was finally selected for several reasons. First, it has been widely used and could be applied to most historical and current satellite sensors, thus allowing consistent long-term studies.
Development and evaluation of seagrass algorithms Algorithm development and evaluation was performed over a large range of NDVI cores (0. 12-0. 77), SPC (0-100%), and SB (0-175. 18 g DW m -2 ) using the Bourgneuf Bay merged dataset (Fig. 3 ). The dataset was randomly split into two groups: one for algorithm development (60% of data), and one for algorithm evaluation (the remaining 40%). SPC data from the two other sites were also used for algorithm evaluation, and its performance was then evaluated for the three sites independently. While a linear relationship was found between NDVI cores and SPC cores (Fig. 3b ), a non-linear relationship was found between NDVI cores and SB cores (Fig. 3c ). An exponential fit (i. e. , NDVI cores = a • (1 -e (b(SBcores+c)) )) was obtained (R 2 = 0. 92; p < 0.
S2 image acquisition and processing Geolocated Level-2A Sentinel-2A/B images of Bourgneuf Bay were downloaded from the European Space Agency (ESA) data portal. Level-2A data have already been atmospherically-corrected using the Sen2Cor processor algorithm  (Main-Knorn et al. , 2017) , and were distributed as bottom-of-atmosphere reflectance (R sat-S2 , dimensionless). We selected only cloud-and shadow-free low tide images (water level at the harbour of Noirmoutier Island, L'Herbaudière < 3. 20 m of the Lowest Astronomical Tide (LAT)), reducing the final S2 dataset to 42 images (see Table  A2  for details). All satellite data processing was performed using the Sentinel Application Platform (SNAP; http://step. esa. int). The atmospheric correction performance was evaluated using two S2 scenes and concomitant in situ reflectance of three types of targets: bare sediment, full seagrass cover, and mixed substrate.
The latter included substrates with bare sediment, seagrass, and/or macroalgae. Different strategies were used for the different target types to obtain the best possible match-ups. For homogeneous areas, namely the areas of bare sediment and dense seagrass cover, the average of the R insitu measurements was compared with the average of several pixels identified as these substrates in the field. For validation over mixed-cover areas, three pixels coinciding with the coordinates of 20 R insitu samples were extracted, and the mean reflectance of those pixels was compared with the mean of the in situ measurements. NDVI S2 was computed for the entire S2 dataset using the bands centred at 665 and 842 nm. To select only the intertidal zone, we applied a geographical mask based on the bathymetric map provided by the French Naval Hydrographic and Oceanographic Service (SHOM).
The lower limit of the selected area was set to 0 m LAT, while the uppermost limit was restricted to 4. 5 m LAT  (Barillé et al. , 2010). Areas of rocky substrate covered by macroalgae near seagrass beds were excluded based on Geographical Information System (GIS) data  (Barillé et al. , 2010). Finally, pixels outside the 0. 12-0. 80 NDVI S2 range were masked. The lower boundary corresponded to either bare sediment (Fig. 3b ) or contamination by a layer of water. The upper boundary was estimated from our radiometric measurements as the saturation of NDVI cores over dense seagrass cover, plus 5%, which was added to create the upper biomass limit mask to account for very dense biomass that can occur in the field, but which was not included in our dataset. According to  Barillé et al. (2010) , values higher than 0.
Statistics for algorithm performance evaluation and map validation Algorithm performance, in situ and satellite product match-ups, and comparison between S2 and WV02 products were evaluated using the regression coefficient of determination (R 2 ), slope of the linear regression, RMSD (Eq. (  3 )), bias (Eq. (  4 )), and mean absolute difference (MAD; Eq. (  5 )), where N is the number of observations, and x corresponds to modelled or observed data. All statistics were computed using the MATLAB software. = = N RMSD (x x ) 1 i N i i 1 modelled, insitu, 2 (3) = = N bias (x x ) i N i i 1 modelled, insitu, (4) = = N MAD |x x | i N i i 1 modelled, insitu, (5).
Evaluation of percent cover and biomass algorithm performance The performance of the SPC algorithm was assessed using the independent in situ datasets collected from three regions along the European Atlantic coast (Bourgneuf, Marennes-Oléron, and Cadiz Bays). SPC cores was predicted from NDVI cores with very good accuracy (RMSD < 5%, R 2 ≥ 0. 98; Table  2 ). Performance was consistent across all sites, suggesting that the same algorithm can be used to seamlessly map seagrass cover over a large geographic range, from 36°N (Cadiz Bay) to 47°N (Bourgneuf Bay). Quite good accuracy was also obtained for SB prediction (RMSD = 5. 31 g DW m -2 , R 2 = 0. 88; Table  2 ), despite the limited number (N = 12) of samples available for evaluation corresponding to SB cores ≤ 51. 4 g DW m -2.
S2 characterization of the Z. noltei seasonal cycle In Bourgneuf Bay, the NDVI S2 exhibited a consistent seasonal cycle characterized by a late summer maximum and a winter minimum (Fig. 7 ). Such a seasonal cycle was expected for this seagrass species in the temperate North Atlantic  (Vermaat and Verhagen, 1996). The dataset combined S2 images from two different orbital cycles (#137 and #94) and two sensors (onboard S2A and S2B). We did not observe any significant difference between the orbits or between the sensors, suggesting that the S2 constellation is consistent enough to be combined to characterize the seasonal dynamics of intertidal seagrass meadows. We then applied a Gaussian model to the NDVI S2 time-series to better appraise the phenological cycle of Z. noltei in Bourgneuf Bay (dashed lines in Fig.
However, due to the limited number of historical images per year, the seasonal cycle cannot be described with the same temporal resolution as we have done with S2 in this work (see paragraph 4. 3). Despite this possible source of variability, an interesting perspective would be to construct a 30-year timeseries (1990-2020) starting with SPOT imagery  (Barillé et al. , 2010)  and continuing with S2 (this study), all the more as this kind of continuity was expected from the launch of the S2 mission  (Hagolle et al. , 2015). Such a spatially-rich and long-term time-series could represent a valuable dataset for environmental management programs, such as the European WFD  (Papathanasopoulou et al. , 2019).
These parameterizations were obtained from our hyperspectral in situ dataset by simulating NDVI values for different sensors using their respective SRF  (Cundill et al. , 2015; González-Audícana et al. , 2006; USGS, 2018). For consistency with a previous study based on the long-term analysis of SPOT imagery  (Barillé et al. , 2010) , we redefined the NDVI thresholds to distinguish two classes of seagrass cover: sparse cover (20% < SPC < 50%) and dense cover (SPC ≥ 50%), and provided consistent thresholds for all sensors (Table  4 ). To complement the NDVI inter-calibration, the difference in the acquisition geometry of each sensor has to be taken into account. For example, the differences between the field-of-view of S2 and L8 (20. 6° and 15° respectively) need to be compensated through BRDF modelling  (Claverie et al.
e. , ± 15 days), and because of the influence of water height on the interpretation of low tide images, our results call for caution in satellite image selection for this intertidal habitat, providing instructions to perform unbiased studies. The S2 dataset is quite recent, but could be complemented by multispectral satellite time-series to investigate long-term changes in seagrass dynamics. For this purpose, we provide guidelines to intercalibrate a multi-sensor NDVI database, and recommend the application of consistent atmospheric correction, if possible, so as to avoid instrumental biases and misinterpretation of temporal changes.
First, the 0% and 100% references were obtained from the average of > 3 measurements over the pure substrates, corresponding to bare sediment and full seagrass cover respectively. Then, the intermediate reference spectra were computed through a linear combination of the references of these two pure substrates. Second, the reflectance dataset was clustered into 11 classes by seagrass cover (0-5%, 5-15%, …, 85-95%, 95-100%), and the corresponding reference spectrum was used to apply the MSC correction to all in situ spectra within each cluster. The radiometric MSC correction applied to in situ spectra clearly improved the quality of the data (Fig. A1a, b ). Also note that, even though this correction was applied independently to each dataset, it improves the significance of ANCOVA tests performed between datasets collected in Bourgneuf Bay in 2018 and 2019.
Acknowledgments We acknowledge the students of the " Ecosystem and Marine Bioproduction " and " Aquaculture, Environment and Society" Master programs of Nantes University  ( 2018/2019  and  2019/2020 ) for their participation in field data acquisition. The  European Space Agency  is acknowledged for the provision of Sentinel-2 and Worldview-2 images, and for the distribution of the SNAP software. We thank the  Laboratory of Planetology and Geodynamics (LPG) of Nantes University  for the use of the ASD spectroradiometer. We also thank the anonymous reviewers for their valuable contributions to improve our manuscript. We especially thank  Dr. Stephanie Palmer  for revising the English of the manuscript.
, 2007; National Ice Center, 2008). Computing lake area change as a function of the distance to the nearest channel required first extracting the Delta Channel Network (DCN) and second computing the distance of every pixel (land or water) to the nearest channel. Automatic channel network extraction for river deltas has recently been advanced with the Python package RivGraph  (Schwenk et al. , 2020) , which utilizes water coverage imagery to extract and skeletonize the DCN. We found that no major channels avulsed or lake breaching took place on either delta during the period of record and therefore used a constant DCN (see the supporting information for more details). To reduce the effect of channelized flow, only water bodies (e. g. , inundated depressions, ponds, and lakes, hereafter collectively referred to as lakes) disconnected from the channel network were considered in our analysis.
This object-based approach, however, can only be performed on imagery with sufficient data quality (i. e. , with a negligible fraction of no data pixels). To implement an object-based approach we used an iterative morphological erosion operation (see Figure  2  inset;  Haralick & Shapiro, 1991)  to classify water pixels according to their position in the different lakes. Specifically, we classified them as belonging to an internal perimeter (IP), with IP 1 indicating the shoreline perimeter of a lake, IP 2 the next internal perimeter, etc. We then estimated a monthly lake shoreline shrinkage rate, S IPi , computed as the fraction of water area loss of water pixels in IP i , as a function of d nc. We note that the weighted average of S IPi for all i, with weights proportional to the number of IP i pixels over the total number of water pixels, will converge to S p when all IPs are considered, i.
Acknowledgments L. V. was supported under the  NASA Earth and Space Science Fellowship Program  grant  80NSSC18K1409. E. F. G. , A. T. , and L. V. received support from  NSF  Grants  EAR-1811909 ,  ECCS-1839441 , and the  UK Research and Innovation Global Challenges Research Fund Living Deltas Hub  Grant  NES0089261. J. S. , A. P. , and J. R. received support from  DOE 's  HiLAT  project. We thank  I. Overeem  for discussions at the early stages of this study and  N. J. Pastick  for his constructive comments and suggestions on the manuscript and for some field survey photos of the Yukon delta he provided that corroborated our interpretation of the NDVI patterns. We also thank  C. Guilloteau  and  A. Longjas  for their insight on remote sensing, change detection, and depositional environments.