dataesr/scientific-paragraphs-categorization huggingface dataset:
Using an APMS, large datasets can be acquired over long time series with little disturbance to the colony. In this study, we refine the CPUE index for ecological studies promoted by  Cairns (1987)  and introduced by  Wilson (1992)  and  Grémillet (1997). We have integrated diet composition estimates measured using faecal DNA metabarcoding with foraging trip duration and body mass data from the APMS. We use this method to determine the availability of key prey fish species of little penguins over time and to estimate the biomass consumed by breeding little penguins.
| Foraging success and trip duration All penguins within the study site were equipped with passive identification transponders (PIT tags: Penguin Parade-Allflex, Australia and Radio Tracking Bay-Trovan, United Kingdom), which are injected subcutaneously between the scapulae  (Chiaradia & Kerry, 1999 ). An APMS with weighbridge, located on the main path in and out of each breeding site, reads these PIT tags, recording individual identity, time and mass, thus allowing us to identify and calculate both the foraging trip duration and the body mass change of individuals during each foraging trip. We used the R statistical language and environment (R Core Team, 2013) to extract foraging trip information from the APMS dataset. We sorted the APMS data by PIT tag number and date, then compared each 'out' record with the previous 'in' record to identify distinct foraging trips made by individual birds.
97) and an OTU table was generated using a custom R script. To assign taxonomy to mOTUs, we used different approaches for each marker. For 18S_SSU, we followed the procedure outlined in  Cavallo et al. (2018) , which provided classlevel identification based on matches in the SILVA database  (Quast et al. , 2013). For 16S_Fish, sequences were resolved to the lowest possible level based on matches to sequences of locally distributed species in GenBank using the NCBI Basic Local Alignment Search Tool (BLAST:  Johnson et al. , 2008)  and then we manually curated the output (Table  S1 ). Following taxonomic identification of sequences, samples were filtered to retain only food sequences (i. e. penguin, parasite and contaminant sequences were excluded). Samples with fewer than 100 food sequences were discarded because these either contained insufficient DNA to analyse or were dominated by non-target DNA so that inferences about diet would be unreliable.
| Developing a prey-specific, predator-derived CPUE index (pCPUE) We filtered APMS records so that only birds that returned on the day before each faecal sample collection were considered (Figure  1 ). We calculated foraging success (body mass change in grams) and effort (foraging trip duration in days) per trip. Weighbridge data can occasionally record multiple penguins crossing (recording erroneously high mass) or not detect a PIT tag (giving erroneously long foraging trips). To account for this, we retained only positive body mass change results <420 g (i. e. the maximum change in body mass for a foraging trip between 1 and 14 days long; Salton,  Saraux, Dann, & Chiaradia, 2015)  and filtered out foraging trip durations >7 days, which are rare for breeding little penguins as they can lead to nest desertion  (Chiaradia & Kerry, 1999; Numata et al.
The 16s_Fish primer only amplifies fish DNA and so only measures the fish portion of the diet. Therefore, we needed to multiply the RRA of the 16S_Fish species sequences by the RRA of all 18S_SSU fish sequences (Actinopterygii) in the samples to determine the proportion of each individual fish species within the total diet. The individual species RRA (e. g. sardine) values were scaled by the mean fish RRA for each breeding season (Actinopterygii: 18S_SSU) rather than by monthly means, which could have been affected by secondary predation  (Sheppard et al. , 2005). See calculation example in Figure  1. The pCPUE was calculated only for fish species that accounted for more than 2% of the diet in one or more years. Code to replicate pCPUE and CPUE from weighbridge and diet data is included in the Supporting Information (Appendix S1).
We also included burrow ID, to determine whether repeated measures at this level had any effect on explaining variation in the data. We set the conditional inference tree to split only for p < 0. 05, a conservative measure, since conditional inference trees are generally fit with a splitting structure of p < 0. 1. Then, to explore the relative importance of variables, we built nonparametric random forests using cforest (package party;  Hothorn et al. , 2006b)  with controls (n = 100,000 trees, mtry = 2, replacement = FALSE) set using the cforest_control function. We retrieved the relative importance of variables from this random forest using the varimp function in party with arguments conditional = TRUE, mincriterion = 0. 95 and threshold = 0. 95  (Strobl, Boulesteix, Kneib, Augustin, & Zeileis, 2008).
| Penguin diet The mean RRA of all fish (Actinopterygii) in the diet, identified by the 18S_SSU primer set, was 62% in 2015-2016 and 64% in 2016-2017 (Table  1 ), but individual samples varied widely. Crustacea, Mollusca, Cnidaria and Tunicata were also identified, the mean RRA of which varied over months and years. Fish had the highest RRA over both seasons, but in 2015-2016 there was moderately high RRA of salps, and in 2016-2017 moderately high RRA of calanoid copepods. Using the 16S_Fish primers, 82 individual fish mOTUs were isolated, with ~75 identified to genus or family levels and 50 identified to species level. Of these, 11 fish accounted for ≥2% of reads in the dataset in one or both years (Table  2 ). From September to December 2015, sardine accounted for 70%-86% of the data in each month.
colony consumed an estimated 618. 0 ± 11. 7 (M ± SE) tonnes of fish during the 2015-2016 breeding season, and 682. 2 ± 14. 4 tonnes in the 2016-2017 season. The sardine biomass extracted by little penguins in this colony during the 2015-2016 season was estimated TA B L E 2 Number of samples analysed using the 16S_Fish primer, with the relative read abundance (%) of species that comprised ≥2% of the dataset in each season, split by season and month Overall catch per unit effort (CPUE; a, b) and preyspecific catch per unit effort (pCPUE); (c, d) for the 2015-2016 (.
DATA AVA I L A B I L I T Y S TAT E M E N T Datasets and code to execute the pCPUE method are available in the FigShare repository at the following: 18s_SSU sample library: https://doi. org/10. 26180 /5ea7d 9b786c4e; 16s_Fish sample library: https://doi. org/10. 26180 /5ea7e 0b0295b0; APMS weighbridge data: https://doi. org/10. 26180 /5ea7d a58da29f; Scat collection.
Number of samples analysed by 18S_SSU primer, with the relative read abundance (%) of classes that comprised ≥2% of the dataset in either year, split by season and month 2015-2016 2016-2017 Aug Sep Oct Nov Dec Jan Feb Sep Oct Nov Dec Jan Feb Mar n Samples 8 57 53 49 54 5 28 24 81 67 17 20 9 6 Actinopterygii 79. 5 54. 3 44. 1 62. 2 82. 8 44. 8 a 93. 6 62 54. 8 55. 7 57 49. 9 95. 3 82 Calanoida 0. 2 2. 8 15. 6 6. 7 0. 6 0 0. 8 32. 8 32. 4 25. 6 1. 8 0 0. 4 7. 1 Copepoda (other) 0 0 0 0 0 0 0 0 0 0 34. 3 18. 3 0 0 Harpacticoida 2. 1 7. 1 1. 5 1. 4 0 0 0 0. 1 0. 3 0 0 0. 2 0 0 Malacostraca 0 0 0. 5 12. 4 3. 8 13. 9 0. 3 1. 1 0. 2 1. 5 0. 7 2. 4 0. 1 0. 4 Oikopleuridae 0 0 3 7. 6 7. 7 11. 5 2 0. 8 6. 2 2. 9 2. 5 28. 8 2. 8 6 Salpidae 11. 9 33. 1 24. 4 3. 3 0 0 0. 9 0. 2 0. 5 2. 5 0 0 0 0 ). In October 2015,TA B L E 1a Low Actinopterygii RRA in January 2016 also coincided with high siphonophore RRA, of which sequences were abundant in that month but did not account for ≥2% of the dataset in that year.
(2006; 2008) ;  Hopke et al. (2008). Besides, only a few long-term PM2. 5 and PM10-2. 5 monitoring data are available elsewhere in this region. This work applied the atmospheric model coupling with air pollution model to study the effect of biomass burning and anthropogenic emission on PM10 concentration. We used the Weather Research and Forecasting (WRF) model version 3. 8. 1 to simulate the meteorological conditions in March 2012. The model's initial and boundary conditions are from the Final Analysis Data (FNL) (NOAA/NCEP, 2000). The modeled temperature and wind speed were compared to the dataset from the Pollution Control Department (PCD). To clarify the model capability, statistical metrics such as Index of Agreement (IOA) and Fractional Bias FB were used for the model evaluation. The output from the WRF model was used as meteorological conditions into the HYSPLIT model to find out the long-range transport of regional air pollutants from neighboring countries of Southeast Asia toward northern Thailand.
The similarity theory scheme is used to emulate thermal gradient over the surface responsible for friction velocity and wind over the surface  (Paulson, 1970; Dyer and Hicks, 1970; Webb, 1970; Beljaars, 1940; Zhang and Anthes, 1982). The model spin-up was conducted from 15 -28 February 2012 to reduce the effect from initial conditions. From March 1, 2012 -April 1, 2012, the WRF model was designed to simulate the weather conditions for the HYSPLIT model. The main meteorological variables of the WRF model, i. e. wind (U, V, W), temperature ( T), surface pressure (Psfc) and relative humidity ( RH), were used as input data for the HYSPLIT model.
Results & discussion To quantify long-range transport from Southeast Asia to Thailand, we begin with a model assessment to examine the performance of the model. We then analyzed the hotpot location in South East Asia using the MODIS satellite dataset based on TERRA and AUQA. The next section is to see the impact of mesoscale meteorology on the smog problem in northern Thailand. Finally, we used a backward trajectory analysis in the HYSPLIT model to find the source location of the smog problem in northern Thailand. there are coupled station sites with a PM2. 5 dataset. Here, we compared the modeled results to the PCD dataset at 3 locations, which is based on the complete PM2. 5 dataset as shown in Table  3. We used the statistical indicator, that is. Index of Agreement (IOA) and Fractional Bias (FB) to examine the efficiency of the model.
The linear analysis of the correlations was based on the Pearson correlation (1), while the non-linear correlation was based on the Spearman equation (  2 ). Pearson's correlation coefficient. ( and analogously for , Spearman Equation; (2) di is the difference between the two ranks of dataset, n is the number of datasets. In general, the linear correlation analysis shows a slightly positive wind speed coefficient and a negative temperature coefficient of 0. 082 and (-0. 46) (Table  6 ). Non-linear correlation analysis shows a negative wind speed correlation with (-0. 53) but a slightly positive temperature correlation with 0. 18. Linear and non-linear correlation analyzes indicate that biomass burning is negative in relation to precipitation and carbonaceous aerosol by (-0. 51) and (-0. 67) (Table  6 ).
The WRF model was performed with 1 domain of 20 km grid spacing covers Southeast Asia, some parts of India, and China. Final analysis data was used as a meteorological initial and boundary condition. The model results including temperature and wind speed were compared to the groundbased measurement from the Pollution Control Department (PCD) in northern Thailand. Model's capability is acceptable when compared to observation which was indicated by Index of Agreement (IOA) in ranges of 0. 57 to 0. 79 for temperature and 0. 32 to 0. 54 for wind speed, while the fractional bias of temperature and wind speed were 1. 3 to 2. 5 °C and 1. 2 to 2. 1 m/s. We analyzed the hotspot data from Near Real Time (NRT) Moderate Imaging Spectroradiometer (MODIS) to find out the proportion of hotspot in Thailand and neighboring countries in Southeast Asia.
A similar approach has been previously used for mixed linear models  [15]  and has been recently adapted in SAIGE  [15]  for MLR, but without an evaluation of its capacity to properly estimate SNP effects. The other method, named the Offset method, which bears similarities with the methods of  [11] , consists of first estimating individual effects in a mixed logistic regression model, and then incorporating these effects as an offset in a (non-mixed) logistic regression model. In the sequel, we evaluate the capacity of logistic regression (LR), MLM, and MLR (using either the score test or the methods mentioned above) to properly take into account population structure associated with heterogeneous prevalence. While Chen et al. were interested in geographically distant populations, spanning several countries in South America, part of our work focuses two geographically very close populations in West Africa, using real genotype data from a recent GWAS  [16].
We also use data simulated with a coalescent model  [17] , reproducing the simulations presented in  [14]. We use similar simulations to evaluate the ability of the PQL and of our two methods to properly evaluate SNPs' effects. Additionally, we propose to generalize the categorization of SNPs in QQ-plots presented in  [14]  to variables other than the population strata, including continuous variables as for example the first PC. The interest of this generalization is demonstrated using the same simulations. All methods are implemented in the R package milorGWAS (for mixed logistic regression in GWAS), freely available on the Comprehensive R Archive Network (CRAN).
Computational efficiency of the methods We measured the running times of GMMAT and of the methods AMLE and Offset, implemented in our package milorGWAS, for a sample size n varying from 1000 to 10, 000. The results are reported in Fig. 1. The figure is on log-log scale, showing that the results are in good accordance with the asymptotic complexity of the methods, given in Fig. 1  Running times for sample sizes varying from n = 1000 to n = 10,000. All computations were run on a Intel i7 at 3. 60 GHz with 8 MB of cache memory. Due to the log-log scale, the timing of an algorithm with complexity O(n α ) appears as a straightline with slope α : slope 3 for the initial step, 2 for AMLE and GMMAT, and 1 for the Offset method the Methods section: the time taken by the initial step (fitting the mixed model) growing as n 3 , it appears linear on log-log scale with a slope close to 3.
Similarly, the complexity O(n 2 ) of GMMAT and AMLE result in straight lines with a slope of 2. GMMAT is con- sistently slower than milorGWAS: milorGWAS runs roughly 8 times faster for the initial step, and 6 times faster for the testing step. Part of the difference in the testing step is due to the fact that GMMAT reads data and write results on disk, whereas milorGWAS uses the RAM. The Offset method is in O(n) while the AMLE and GMMAT are in O(n 2 ) , which makes this method appealing in terms of running times. However, in practice the initial step in O(n 3 ) is the limiting factor: it would take more than 24 hours to fit the MLR model for n = 50,000. To manipulate larger samples such as those of biobanks, more complex workarounds and other approximate methods are needed  [15].
Type I error in the presence of population structure We analyzed two data sets with a simulated binary phenotype. Simulations were based either on real genotype data from a recent GWAS  [16] , or on large cohort simulated with a coalescent model. Both data sets present cryptic relatedness and population stratification, with two strata (or two cohorts) with different disease prevalence. The simulations are fully described in the Methods section. The first simulated set uses genotypes of 800 individuals from South Benin, ascertained in two sites distant of 20 km from each other, forming the two population strata in which we set different prevalences for simulating the phenotypes. These data were analyzed using the logistic regression (LR), the mixed logistic regression (MLR) (with Chen's score test GMMAT or equivalently Wald test with AMLE) and the mixed linear model (MLM) (the status being analyzed as a quantitative trait with values 0 or 1), with or without the top 10 principal components (PCs) in the model.
However, both logistic regression and MLR show an adequate correction for all categories of SNPs. The second simulation set consists in 10, 000 individuals simulated on a 20 × 20 grid, whose genotypes were obtained from a coalescent model. These data include some first order relatives. Ten millions independent SNPs were simulated, among which 2,840,903 had a minor allele frequency above 5%. Two strata were defined according to the indi- viduals position on the grid (a "high risk strata" was defined as being the top left quarter of the grid), and a binary phenotype was simulated with different prevalences on these strata. Similar analyzes were performed on these data (Additional file 1: Figure  1 ). The 2,840,903 independent SNPs are at 23, 7% in category 1, 58. 8% in category 2 and 17. 5% in category 3.
20 km apart) from South Benin  [16] , but with different self-reported ethnicities. Heterogeneity of prevalence may result from environmental factors (e. g. lifestyle, nutritional behavior, etc. ), and could occur frequently in association studies, thus making the analysis of binary traits with the MLM incorrect, in particular in populations with a high genetic diversity, such as African populations. A similar result was retrieved with data simulated with a coalescence model on a square grid, the "high-risk strata" consisting in the top left quarter of the grid; these simulations included many first-order related individuals and a random effect based on the kinship matrix. The mixed logistic regression (MLR) including the top PCs as fixed effects is the only method to completely correct for population structure in both simulations.
All methods are available in an R package milorGWAS based on the R package Gaston  [21]  for data manipulation. The R and C++ source code of milorGWAS is available on the CRAN at https ://CRAN. R-proje ct. org/packa ge=milor GWAS. An excerpt of the data simulated using the coalescent model is also available, and is used in the package's vignette to illustrate the methods.
Asymptotic complexity of the methods and efficiency of the implementation The initial step of both AMLE and Offset method is to fit the MLR model logit E(Y ) = Xβ + ω. Each iterative step of the PQL needs the inversion of an n × n matrix, where n is the number of individuals; the complexity of this operation is O(n 3 ). The second step of the AMLE involves, for each SNP, multiplying an n × n matrix P by the vector of genotypes; this is the most costly operation, and the complexity of this step is thus O(n 2 ). The second step of the Offset method is an iterative algorithm, each itera- tion of which is in O(n) (considering the number of covariates as fixed). The complexity of the second step of the Offset is thus O(n). Our package milorGWAS is implemented in C++ using Rcpp  [22] , and RcppEigen  [23]  for matrix arithmetic.
Simulations based on a coalescent model We also performed coalescent simulations, reproducing closely the simulations described in  [14] , to obtain genotypes for a large cohort of 10,000 individuals, with both population structure and relatedness, using the ms software  [17]. This procedure, which is based on a stepping stone model with symmetric migration between adjacent cells of the grid, is commonly used to simulate a population with a spatially continuous population structure  [24, 25]. We use a 20 × 20 grid, in which the migration rate between adjacent cells was set to 10; this parameter produces a Wright's fixation index F st < 0. 01 when dividing the simulated grid into two equal sub-populations, a level comparable to what is observed within Europe  [24]. We simulated a total of 10 million independent SNPs.
After filtering out SNPs with MAF lower than 5% , 2,840,903 SNPs were available. The full command line arguments for ms are included in Additional file 1. We also created a R data package containing part of the simulated data (link in Additional file 1). To obtain related individuals, we first simulated the genotypes for 8000 founders (20 on each of the 400 cells). We then sampled 10 pairs of individuals in each cell, forming 4000 couples, and simulated two offsprings by gene dropping. Thus we obtained 16,000 individuals (founders and offspring) from which 10,000 individuals were randomly selected to obtain the cohort. The phenotype was simulated as before using Eq. (  2 ), where Z i was set to one when individuals were sampled in the top left 10 × 10 grid, corresponding to strata with a higher risk.
Illustration on real data To illustrate the method, we applied it on the data from the GWAS in South Benin described in  [16]. We used as binary phenotype the presence/absence of any malaria infection during the follow-up (there were 229 individuals with no infection and 546 with at least one infection). We tested all SNPs with a minor allele frequency greater than 0. 05 on a 100 kb segment on chromosome 20, which contains one of the strongest association signals discovered in  [16]. We included as covariates the site of ascertainment, the duration of the wollow-up, and mean infection exposure. We performed the testing with AMLE and, for comparison, with the PQL.
Availability of data and materials The genotype datasets from South Benin used in the current study are available in the DataSuds repository (https :// doi. org/10. 23708 /EXSQT M). The dataset are not publicly available due to patient confidentiality but are available for researchers who meet the criteria for access to confidential data. Data to reproduce the simulations based on a coalescent model are included in this published article and its supplementary information files.
Introduction Companies often need to engineer a portfolio of software variants instead of one-size-fits-all solutions. Creating variants allows tailoring systems towards varying stakeholder requirements-different functionalities, but also non-functional requirements, such as performance or power consumption. Variant-rich systems are especially common in traditional engineering domains including automotive, industrial automation, telecommunication. In addition, recent trends such as the Internet of Things (IoT)  [15] , cyber-physical systems  [113, 156] , or robotics  [84] , further increase the need for customization. Software product line engineering (SPLE) aims at effectively engineering a variant-rich system-a software product line-in an application domain. SPLE advocates establishing an integrated software platform from which individual variants can be derived, typically in an automated, configuration-driven process.
We used document analysis, semi-structured interviews, and focus groups on cases that cover a wide range of domains and development scales, from a rather small web-application engineering case to ultra-large software engineering for automotive or industrial component production. In addition, we conducted a lightweight literature review on relevant SPLE adoption case studies, experience reports, surveys, and meta-studies, supporting the formulation and synthesis of challenges we present. Our industrial cases primarily represent the development in a small part of a company, such as a single division or development team. We refer to these as use cases in the remainder. We also investigated cases provided by tool vendors who are looking to integrate product-line engineering concepts into their tools, referred to as tool cases.