dataesr/scientific-paragraphs-categorization huggingface dataset:
Delay Propagation We model the interaction between trains with a mechanism of propagation of delay from other trains to train i. The investigations performed on the datasets suggest that such interaction can occur only when nearby trains are in the Backward propagation configuration of Fig. 3. Following the picture of the model in Fig. 7 , when the train i starts travelling on the link BC leading to the C station is susceptible of propagation from the train j, which is travelling from C towards D. Inspired by the SIS models of epidemic spreading  [24] , in which agents can get infected and recover from the illness with a fixed probability, we introduce the delay propagation parameter β ∈ [0, 1] describing the probability that the delay of the train j propagates to i. We assume that the propagation occurs at the time i starts travelling onto { for t'j = tj 1for t'j = 0 Fig.
For this purpose we provide a definition of "cluster of congested stations" and how to discriminate whether a station is "congested", i. e. when its functioning is inefficient because of the hoarding of delays. For each station we can define a threshold between the "functioning" and the "congested" status based on the value of incoming train delays. We calculated such threshold for each station as the average value of the delays of all the trains arriving at the considered station in the dataset. Considering a period of two months (March and April 2015), we examined each station every 5 minutes during the day, checking whether the average delay of all the trains moving towards that station during that lapse of time was above the average. We consider the station as "congested", while if the average delay is below that threshold we consider it as "functioning".
In this work, we developed a novel model of delay emergence and propagation between trains moving on Railway Networks, partly based on empirical laws inferred from the data governing the accidental emergence of delays from the influence of adverse occurrences. Remarkably, the statistical models used for the description of this "exogenous delays" showed how these adverse occurrences are the result of the same finite number of causes (e. g. , bad weather conditions and malfunctions) independently from the topology, as opposed to the scale of these disruptions, which are largely influenced by the part of the network the trains are travelling on. We find, in fact, that both the complexity of a station as measured by the number of different routes originating from it and the length of a route, are connected to the magnitude of the microscopic delays composing the overall delay by simple statistical laws with two parameters, whose value is the same all over the network and can be inferred and fixed by data.
Moreover, the modelization of exogenous delays could be refined by taking into account more static topological feature of the Railway Networks or by introducing dynamical variations due to different traffic conditions. Another interesting possibility would be to increase the number of Railway Systems under study (e. g. , the French SNFC sharing similar characteristics with respect to the Italian and German systems whose data is publicly available  [36] ), in order to understand the reasons behind the quantitative differences observed, e. g. larger delays in the Italian case. In the final part of the paper, we show how the model can be applied to study the capability of functioning of the system after a localized large disruption occurring in a single node of the network. This approach can be easily extended to more complex case studies of distributed disruptions, the occurrence of strikes, the removal of trains or parts of the network, also including the possibility to introduce delay management strategies to increase the resilience of the system.
Hence, we downloaded all the information displayed on the website each day at 23:30 in order to be sure that each train would have arrived at destination. Starting from the 1st of January 2015 and for the whole 2015, we collected 12 months of historical data about the dynamics of regular and high-speed trains in Italy. For each train we get an identifier, the ordered list of stations the train has to cross, the scheduled arrival time at each station and its delay. The resulting dataset comprehends the traffic running on 2253 stations, with a daily average schedule pertaining 8112 trains on 7062 links. Note that "ViaggiaTreno" does not collect information about the geographical position of the stations. Such information has been integrated by means of Wikipedia and Google Maps, allowing us to represent the geo-localized network of Italian railways.
A. 2 The German Dataset The data about German Railways have been collected through the OpenDataCity 5  website. This site gathers different datasets collected by a variety of on-going or terminated projects dealing with open data. In particular, the data we analysed come from the "Zugmonitor" project, which aimed at providing a web-app and an API to German travellers in order to have real time information on the position of the trains on the German Railway Network and their delay. The project is no-longer running and the API is not accessible anymore. However, some dumps of historical data collected during the project are still available. In particular, we downloaded all the data regarding year 2015, covering the same period of the Italian dataset. This dumps collected not only the delay at each station like in the Italian case, but also the delay at intermediate points between two stations.
All the points are also geo-localized so that it is possible to reconstruct a quite accurate trajectory of the trains. In order to be consistent with the Italian dataset we used the geo-localization only to identify the position of the stations in the map. Scheduled arrival times at each station were also stored in the dump, so that in the end we managed to reconstruct a dataset with a structure identical to the Italian one. The resulting dataset includes data for 5979 stations with a daily average schedule containing 11,975 trains on 16,277.
B High-speed layer The structure Italian and the German Railway Networks is the overlap of two distinct layers, the normal-speed and the high-speed one. These two layers are different from the structural point of view. The high-speed layer in fact has to allow for fast travelling trains and have a different kind of rails connecting stations and, in general, it is reasonable to assume that highspeed trains and normal-speed ones do not interact when travelling from a station to another. However, the nodes of the network (i. e. stations) are shared between the layers, making the network a "multiplex"  [40, 41]  and allowing for interactions of the two different kind of trains. Our datasets contain information about high-speed and normal-speed trains, for both the Italian and German case and in principle it could be possible to study the dynamics of the high-speed layer and its interaction with the normal one.
In the main text we decided though to focus on the normal layer cutting-out the high-speed part. This choice was made for sake of simplicity, since the rules of interaction between the layers might have been hard to understand or derive with data analysis. Here, we will show that this approximation is reasonable due to the smaller numbers of the high-speed trains travelling and their poor effects on the dynamics of the normal-speed one. In our datasets it is possible to identify high-speed trains thanks a specific identifier ("ES*" for the Italian Network and "EC", "IC", "ICE" for the German Network) and use them to build the High-speed layer of the Railway Network in a similar way that has been done for the normal layer in the main text. The number of travelling high-speed trains per day is considerably smaller with respect to the normal-speed trains, being of ∼ 210 and ∼ 1055 for the Italian and German case respectively.
D Departure Delay Distributions An approach similar to the one adopted for the exogenous delays on the links of the networks can also be adopted for the departure delays at the stations. In this case we categorized the departing stations (i. e. a subset of the nodes in the network) according to their out-degree kout whose distributions are shown in Fig. K. Having divided the nodes of the network according to kout, we can fit the aggregated departure delay distributions as kout varies as shown in Fig. M , L and N. Note that negative departure delays are not present in the German dataset. These distribution have been fitted using a q-exponential functional form as in equation 4. The behavior of the q and b parameters for these distributions according to kout can be summarized by the equations: q(kout) = conts, b(kout) = Ae -ak out.
(6) The parameters for equations (  6 ) are obtained by fitting the empirical data as shown in the main text. 6 ), governing the behavior of the parameters q and b of the q-exponential distribution as out-degree kout of the nodes varies. Similarly to what we have found for the dependency of the exogenous delay with the length of the links, in the German Railway Network no dependence can be observed and the probabilities of having a positive or negative departure delay can be considered constant in every station. However, this is not true for the Italian Railway Network where just the probability of having a negative delay can be considered constant. On the contrary, the probability of having a positive departure delay increases linearly with kout.
E Optimal Choice of β In order to determine the optimal value of the β, we tune our model to reproduce with the highest probability the delay that each train gets whenever it crosses a station during its path. Considering a train i arriving at a station n on a given day, we call δt emp i,n its measured arrival delay at that station as recorded in the dataset. Hence, we perform 200 simulations of the schedule of the considered day in order to compute the distribution P (δt i,n ) of the corresponding δt i,n. Hence, considering the null hypothesis that δt emp i,n is produced by our model (i. e. it is extracted by P (δt i,n )), we calculate the double tailed p-value for such hypothesis for every pair (i, n), i. e. for every train and for every station. For each day we average the p-values of all the train-station pairs to obtain the performance metrics for β.
Fig. R  shows the distribution of the average weight of the links (in terms of the number of trains that have travelled in the links in our whole datasets) and the distribution of the length of the links connected to a station, in the case of stations with p-value larger and smaller than 0. 05. We can see that in the latter, the distribution of the weight is considerably narrower, indicating that in these stations the traffic is usually low. Hence, the disagreement might be the result of a poor sampling of the exogenous disturbances around these stations leading to poor predictions or to a dependence of the transmission parameter β on the traffic conditions that have been ruled out when we assumed them to be constant in time and uniform all over the network. For the Italian case, these stations also are usually connected to links with a shorter distance.
Funding This work has been supported by the  KREYON  Project, funded by the  John Templeton Foundation  under contract n. 51663; VDPS acknowledges financial support from the  Austrian Research Promotion Agency FFG  under grant # 857136. RDC as Newton International Fellow  of the Royal Society acknowledges support from  The Royal Society ,  The British Academy  and the  Academy of Medical Sciences  ( Newton International Fellowship ,  NF170505 ).
However, the 3D scan and µCT scan datasets performed identically in classifying individuals based on intra-specific patterns of sexual dimorphism. Discussion. Compared to µCT scans, we find that even low resolution 3D scans of very small specimens are sufficiently accurate to classify intra-specific differences. We also make three recommendations for best use of low resolution data. First, we recommend that extreme caution should be taken when analyzing the asymmetric component of shape variation. Second, using 3D scans generates more random error due to increased How to cite this article Marcy et al. (  2018 ), Low resolution scans can provide a sufficiently accurate, cost-and time-effective alternative to high resolution scans for 3D shape analyses. PeerJ 6:e5032; DOI 10. 7717/peerj. 5032 landmarking difficulty, therefore users should be conservative in landmark choice and avoid multiple operators.
Photogrammetry, like 3D surface scanning, is a low-cost alternative to µCT and comes with its own trade-offs in time and scan resolution  (Katz & Friess, 2014). Compared to the new generation of blue light surface scanners, photogrammetry requires more time for image acquisition and for file processing  (Katz & Friess, 2014). A previous study on a single macaque specimen reported inconsistent levels of error across operators and scanners, which contributed to the lack of general pattern for differences across scanners/resolutions  (Shearer et al. , 2017). However, using an interspecific dataset,  (Fruciano et al. , 2017)  reported higher repeatability for the higher resolution scans and 2. 07-11. 26% of total variance due to scan type (depending on device, operator and landmark set combination).
Data collection We selected 19 adult individuals, male and female, of Pseudomys delicatulus from the Queensland Museum in Brisbane, Australia (specimen numbers and sexes in Table  S1 ). The cranium from each individual was scanned at the Centre for Advanced Imaging at the   How does variation due to scan device compare to other sources? Table  1. and when variation due to bilateral asymmetry is removed? Table  2  How does variation for the symmetric shape component look? Does repeatability (i. e. operator error) di er by scan device? Table  3  Does repeatability di er if fewer landmark types are used? Table  3  D.
The scanner was operated at 80 kV energy, 250 µA intensity with 540 projections per 360 • , a medium-high magnification with bin 2 was applied, and 2,000 ms exposure time. The samples were scanned at a nominal isotropic resolution of 28 µm. The data were reconstructed using a Feldkamp conebeam back-projection algorithm provided by an Inveon Acquisition workstation from Siemens (IAW version 2. 1, Munich, Germany). Surface models were obtained using Mimics Research version 20. 0. Each cranium was also scanned by a HDI109 blue light surface scanner (LMI Technologies Inc. , Vancouver, Canada) with a resolution of 80 µm. For brevity, we will refer to this method as 3D scanning. For this method, the cranium was placed on a rotary table providing the scanner with 360 • views. To capture the entire shape, the cranium was scanned in three different orientations: one ventral view with the cranium resting on the frontals and two dorsal views with the cranium tipped to each side, resting on an incisor, auditory bulla, and zygomatic arch.
To assist others in replicating our HDI109 3D surface scanning on small specimens, we have included a Standard Operating Procedure with our settings (Supplemental Information 1). After scanning every individual with both scan methods, we then replicated each 3D model three times so that each individual was represented by six replicates, giving a total sample of 114 3D models (Fig. 2A ). Each 3D model was landmarked in Viewbox version 4. 0 (dHAL software, Kifissia, Greece; http://www. dhal. com; Polychronis et al. , 2013). To capture shape, we placed 58 fixed landmarks, 145 semi-landmarks on curves, and 86 patch points (points that during sliding are allowed to slide across a 3D surface defined by the 3D model and semi-landmark borders) for a total of 289 points (Fig. 3 , Table  S2 ). We used the template feature in Viewbox to semi-automate the placement of semi-landmarks on curves and to fully automate the placement of patch points.
Data analysis The landmark coordinates for all 114 3D models were aligned using a generalized Procrustes analysis followed by projection to the tangent space, as implemented in the R package geomorph (v. 3. 0. 5) (Adams, Collyer & Sherratt, 2016; Adams & Otarola-Castillo, 2013). Generalized Procrustes analysis of each set of landmark coordinates removes differences in size, position, and orientation, leaving only shape variation  (Rohlf & Slice, 1990). Semi-landmarks and patches were permitted to slide along their tangent directions to minimize Procrustes distance between 3D models  (Gunz, Mitteroecker & Bookstein, 2005). The resulting Procrustes tangent coordinates were used as shape variables in all subsequent shape analyses. All our statistical analyses were performed either in R (v. 3.
This value is normally between 0 and 1, with values close to 1 indicating much larger variation due to the factor used in computing the Procrustes ANOVA (in our case, variation among individuals) compared to residual variation (in our case, variation among digitizations). In other words, comparing repeatability between devices gives similar information to that obtained by the box plots of Procrustes variance but on a more easily interpretable scale from 0 to 1. We repeated our computations of repeatability for subsets of the data to test whether introducing semi-landmarks on curves and surfaces (patch points) changed the repeatability relative to a fixed landmark-only dataset. We did so for both 3D and µCT datasets to see if these trends differed by scan device. Finally, we investigated whether there is a difference between devices in a common task: the correct classification of sexual dimorphism (Fig.
For the 3D dataset, the error term (Res/Rep) contributes 17. 8% of variation while asymmetry (Side) contributes 19. 2%. In other words, for our 3D scan dataset, error contributes almost as much variation as asymmetry (Table  1B ). The Procrustes ANOVA for just the µCT dataset, however, did not have this problem to the same degree. Here, the error term (Res/Rep) contributes only 8. 52% of variation while asymmetry (Side) contributes 20. 8% (Table  1C ). In other words, error contributes less than one half of the contribution of asymmetry in the µCT dataset. The Procrustes ANOVA on just the symmetric component of shape reports the individual shape variation, representing biological variation, is 73. 4% (Table  2 ). Differences between scan devices represent 14. 6% and the residuals encompassing differences among replicates or operator error represent 12.
Procrustes variance and repeatability To compare the digitization error in each scanning device dataset, we calculated the Procrustes variance among the replicate triads of each individual. We found that Procrustes variance is significantly (p < 0. 001) higher in 3D scans (mean = 1. 31 ×10 -4 ) than in µCT scans (mean = 4. 76 ×10 -5 ) (Fig. 6 ). This means that digitizations are more variable in 3D scans than in µCT which is consistent with decreased clustering in 3D scans relative to µCT scans in the PCAs (Fig. 4 ). The repeatability for each scan dataset mirrored the Procrustes variance results. We found that the µCT scan dataset had a repeatability of 0. 896 and the 3D scan data had a repeatability of 0. 750 (Table  3A , Table  3D ). This means operators are more successful at repeating their digitizations (i.
e. , landmark placements) with µCT scans than with 3D scans. To test how different types of landmarks impacted repeatability, we calculated repeatability for combinations of landmark types for 3D and µCT datasets consisting of only the symmetric component of shape (Table  3 ). Because sliding landmarks depend on the placement of fixed landmarks (and patch points depend on both fixed and semilandlandmark curves), we could not isolate each type of landmark's repeatability. The Table  3  Comparison of operator error in 3D scan and µCT scan datasets using Procrustes ANOVAs and repeatability scores. The repeatability score (R) is a value that reflects the ease of digitizing in a repeated measure study design. It is calculated from the Procrustes ANOVA using formulas for the intraclass correlation coefficient.
The Procrustes ANOVAs were found by subsetting the dataset by scan device and by landmark types and then performing separate generalized Procrustes and bilateral symmetry align-  analyses restricted to completely manually placed fixed landmarks always had the lowest repeatability of the three types of landmarks (Table  3C , Table  3F ). Repeatability was always highest for the datasets including all three types of landmarks including the semi-automated semi-landmarks and the completely automated patch points (Table  3A , Table  3D ). Higher repeatability in datasets with the sliding landmarks may result because the sliding smooths out user placement error across replicates.
Analyses with a biological example: sexual dimorphism A small subset of our dataset had sex information (n = 11; f = 7, m = 4), allowing us to perform a test of whether using different scan devices classify males and females according to shape with the same level of accuracy. Our Procrustes ANOVA on the symmetric component of shape variation using sex and device as factors found that shape differences due to device (Rsq = 0. 0646) and sex (Rsq = 0. 0952) are both significant (p < 0. 001). Both factors have relatively small effect sizes, however, sex captures slightly more shape variation than device (Table  4 ). However, the between-group PCAs do not suggest marked sexual dimorphism to begin with (Fig. 7 ). Therefore, the subtlety of this biological signal could be the main reason for the small contribution of sex to total variation.
DISCUSSION In this study, we contrasted very high resolution µCT scans with their extreme opposite: 3D surface scans of very small specimens. Our low versus high resolution datasets allowed us to assess whether the low resolution scans still allow defensible investigations of biological shape variation. We found that despite the low quality of the 3D scans, sufficient amounts of biological variation are present to perform, at the very least, typical interspecific comparisons. In datasets with only very slight intra-specific differences, more difficulties in distinguishing biological signal from the noise introduced by error during data collection. For example, the subtle sexual dimorphism in our small sample was only just distinguished. However, we present three considerations to make before using low resolution datasets.
First, we found that variation due to scan device and digitizations is substantial relative to asymmetric variation. This makes low resolution datasets a poor choice for studies on asymmetry. Second, using 3D scans creates more random error due to increased Each PCA provides an exploratory visualization of shape variation between males and females in our subsample with sex information (n = 11). Males (n = 4) are plotted in light blue and females (n = 7) are plotted in dark red. Results from the cross-validation test can be found in Table  5. Full-size DOI: 10. 7717/peerj. 5032/fig- 7  Table  5  Between group PCA classification test to assess whether one scan device dataset performs better at identifying sexes based on shape. This analysis averages shape among replicates, computes a between-group PCA separately for 3D and µCT datasets, and runs a cross-validation classification test.
The results indicate whether one type of scan dataset is more successful at classifying males versus females based on the shape variation present in the dataset. It also returns a kappa statistic; a kappa value over 0. 20 indicates ''fair'' agreement between the two datasets. Shape variation visualized by sex can be seen in Fig. 7. landmarking difficulty, therefore care should be taken in landmark choice, and possibly landmarking software and operator choice. Digitization error may also be reduced by taking averages of repeated measurements  (Arnqvist & Martensson, 1998; Fruciano, 2016). Third, using 3D scans also introduces a source of systematic error relative to µCT scans, therefore we recommend not combining them whenever possible (see also  Fruciano et al. , 2017) , and especially in studies on small intra-specific variation.