dataesr/scientific-paragraphs-categorization huggingface dataset:
com/auranic/ClinTrajan and in the Section 'Method of Elastic Principal Graphs (ElPiGraph)' of this manuscript. The principal tree explained. % of total variance in contrast to the first two principal components that explained. % and first five PCs explaining %. The obtained principal tree (shown in Figure ) was used to provide a D layout of the dataset which can be used for visualization of various clinical variables and the results of analyses. Globally, the principal tree defined three terminal non-branching segments populated with non-lethal clinical cases (indicated as # ,# ,# in  ingly. Each node of the principal tree is connected with a subset of data points. We performed the enrichment analysis, based on application of independence χ test in order to determine the node which is the most strongly associated to  'no complication' class (black points in Figure ).
Assigning data point classes As we've mentioned above, the classes of the clinical observations can be usually defined by selecting a subset of clinical variables which represent some final read-outs of a patient state. Thus, in the myocardial infarction complications dataset, clinical variables report the complications, of them represent binary variables and categorical variable LET_IS, whose value is if there is no lethal outcome. Otherwise, LET_IS can take one of the nominal values representing the death cause. Following the methodology suggested in this study, the LET_IS variable is first made a subject of dummy coding, introducing binary features. The resulting binary variables were characterized by unique combinations of / values, which looked too many to define one class per such unique combination.
Non-lethal class ' ' was distinguished by 'P_IM_STEN' (post-infarction angina), and the class ' ' -by the presence of diagnosed 'A_V_BLOCK' (third-degree atrioventricular block). Besides this clustering, a particular non-lethal state was distinguished characterized by zero values of all complication variables. We distinguished this class as a separate 'no complications class'. In the complete dataset, it corresponded to % of clinical records (denoted as black points in Figure  and Figure  ). In the rest of the analysis, all complication variables have been analyzed together with the clinical characteristics.
Dataset partitioning (clustering) by principal tree non-branching segments The explicitly defined structure of the computed elastic principal tree allows partitioning the dataset accordingly to the projection of the data points on various internal and terminal segments as described in the 'Methods' section and shown by color in Figure , panel 'Tree segments (branches)'. Such partitioning can play a role of clustering with the advantage of that the tree segments can recover non-spherical and non-linear data clusters. In addition, the data clusters, defined in such a way, are connected in a treelike configuration, with junctions corresponding to the branching points which can correspond to 'points of no return' in the state of the patients. Each non-branching segment in the tree can be associated by.
If a data point (clinical observation) is attributed to several trajectories then it is characterized by the same pseudo-time value on each of them. This can be interpreted as the state of uncertainty from which several clinical scenarios can be developed in the further course of the disease, following one or several bifurcation points. Those clinical observations belonging to a single trajectory correspond to less uncertainty in the prognosis, with higher chances to end up in a terminal state. In order to determine the factors affecting the choices between alternative clinical trajectories, it is necessary to associate clinical variables to each of the trajectory and determine the trend of their changes along them. Mathematically this corresponds to the solving the regression problem connecting a clinical variable and the observation pseudo-time.
Predicting survival and lethal risk factors Each clinical trajectory extracted from the analysis of synchronic clinical data is interpreted as a possible ordered sequence of states from the least heavy condition to the extreme final point of the trajectory. Assuming that for a given patient state all the downstream points on the clinical trajectory represent possible future states of the patient, we can make a prediction of possible clinical risks connected with moving along this trajectory. In particular, this can be used for estimating lethal risks, is such events are recorded in the clinical dataset. In order to evaluate the risks of a clinical event in the future, well-developed methodology of survival analysis can be used, but using the pseudotime value instead of the real time value.
We will call such analysis the pseudotime survival analysis. The pseudotime quantified along different clinical trajectories might be incomparable in terms of the physical time. Therefore, the pseudotime survival analysis should be performed for each trajectory individually, even if the estimated risks can be visualized together using the common pseudotime axis. We applied a non-parametric estimator of the cumulative hazard rate function (see Methods) in order to quantify lethal risks along ten identified clinical trajectories in the myocardial infarction complication dataset (Figure  , A ). This analysis highlighted out of trajectories as characterized by elevated hazard rates of lethality, which is a quantification of the distribution of lethal cases on the principal tree shown before in Figure.
Clinical trajectories in large-scale observational diabetes data In order to check if the ClinTrajan package can be applied to larger datasets, we extracted clinical trajectories using a publicly available dataset, representing years ( -) of clinical care at US hospitals and integrated delivery networks. The dataset contains records satisfying the following conditions: ( ) it is an inpatient encounter (a hospital admission), ( ) it is a diabetic encounter, that is, one during which any kind of diabetes was en-tered to the system as a diagnosis, ( ) the length of stay was at least day and at most days, ( ) laboratory tests were performed during the encounter, ( ) medications were administered during the encounter. The data contains such attributes as patient race, gender, age, admission type, time in hospital, medical specialty of admitting physician, number of lab test performed, HbA c test result, diagnosis, number of medication, diabetic medications, number of outpatient, inpatient, and emergency visits in the year before the hospitalization.
In the supervised setting, the aim of the analysis of this dataset is usually to predict the readmission event ('readmitted' variable) within days after disposal from the hospital. In our analysis, we considered the readmitted variable as a part of the data space, in order to perform unsupervised analysis of the dataset with the aim to extract clinical trajectories, some of them leading to the increased readmission likelihood. The exact protocol for encoding the diabetes dataset is provided at the ClinTrajan github https://github. com/auranic/ ClinTrajan/. Importantly, we encoded several categorical variables as ordinal. In particular, the readmitted variable was encoded in three levels with value corresponding to 'No' (absence of recorded readmission), -to '> days' and -to '< days'. The 'A Cresult' feature (related to the HbA c test) was encoded in two variables.
% of records corresponding to the elapsed states (hence, without possibility of readmission) have been excluded from the analysis, similarly to some previous analyses [ ]. The resulting encoded dataset contained variables ( numerical, ordinal and binary). We performed the data pre-processing similar to the one done in the infarction datasets, with imputing the missing values in the 'A Cresult_value' column and with further application of optimal scaling to ordinal values (see the corresponding Jupyter notebook at https://github. com/auranic/ ClinTrajan/). The dimensionality of the dataset was reduced to as it was the consensus value resulted from application of several methods of intrinsic dimension estimation (see Supplementary Figure ), excluding outlying measurements. The principal tree algorithm was applied with the same parameters as in the previous section.
The construction of the principal tree with nodes for the -dimensional dataset with data points took approximately seconds on an ordinary laptop. The principal tree explained % of the total variance in contrast to % of the first two principal components, with PCs needed to explain the same percentage of variance as the principal tree. The tree contained branching points (see Figure  , A ) with one forthorder star. The principal tree-based data layout was used to visualize the values of data space variables (Figure  , A ), and some other variables from the annotation data (Figure  , B, D ), some of which did not participate in determining the structure of the principal tree. As a root node in this case, we selected the middle node of one of the internal segments of the principal tree (segment # in Fig-  ure , A ), which was characterized by the shortest times spent in the hospital, smallest number of all procedures, no history of inpatient stays or emergency calls in the preceding year, normal predicted (not measured) value of HbA C, absence of any medication.
Trajectory-based analysis of the relation between early readmission rate and the measured glycated hemoglobin HbA c In the original publication of the diabetis dataset, several observations have been reported [ ]. First, it was observed that the dependency of early readmission (in less than days) frequency estimate on the fact of the measurement of HbA c is conditional on the type of the primary diagnosis (with three major ones being di-abetis, circulatory and respiratory diseases). Second, for the patients with diabetis as primary diagnosis, it was shown that not measuring the level of HbA c is connected to the increased risk of early readmission. Interestingly, from the Figure  of [ ],  one can conclude that, in patients with diabetis as primary diagnosis, high levels of measured HbA c are connected to decreased readmission risk compared to the normal level of measured HbA c.
This quite paradoxical observation was done in simple calculations of the early readmission frequency as well as in rate calculations adjusted for several clinical covariates. In order to illustrate the advantage of the trajectory-based patient stratification, we recomputed the simple not adjusted estimations of the early readmission frequency as a function of measured HbA c in sets of patients with different primary diagnosis (Fig-  ure , A ). This reproduced the previously made conclusions from the original study [ ]. The frequency of readmission appeared to be higher in the patients with not measured HbA c. Qualitatively similar to the previous publication, the readmission rate was significantly lower for the high values of HbA c compared to its normal levels (. % vs. %). Note that the analyzed dataset has changed since it original publication with more than new patients being added.
Discussion In this study we considered two rich and large publicly available observational clinical datasets from the two most challenging areas of public health: cardiology and diabetes. Both datasets contain syncronic (related to the moment of staying in the hospital) observations over a relatively large population of patients. Therefore, the traditional unsupervised machine learning approach for treating these data in order to classify clinical states is supposed to be some kind of clustering or manifold learning. We demonstrate that there exists an alternative approach which allowed us to represent these data as pseudo-diachronic, i. e. , reflect to some extent the temporal aspects. This opens a possibility to classify not only the states of particular patients but their hypothetical clinical trajectories arriving from the past and projected into the future.
This in turn gives a possibility to reason in terms of dynamical disease phenotyping, e. g. , classifying clinical states in terms of the disease dynamics type. Identification of clinical trajectories is made possible due to the use of the branching pseudotime approach, consisting in modeling the geometry of the dataset as "bouquet" of diverging trajectories, starting from one or several hypothetical quasi-normal, e. g. , characterized by the least severe condition, disease states. The progression along a particular clinical trajectory can be quantified in terms of pseudotime, reflecting the abstract accumulated amount of changes in the observed clinical traits. The main requirement for the possibility of such reconstruction is the existence of sufficient number of observations (thousands) such that the individual variations in the clinical states would reveal the major non-linear routes along which they progress in real physical time.
This effort resulted in ClinTrajan Python package which can be readily applied in the analysis of clinical datasets containing even millions of observations. In the real life diabetes dataset considered here and containing more than one hundred thousands of observations, the analysis by ClinTrajan takes few minutes on an ordinary laptop. Use of ElPiGraph is the most relevant in the case when the hypothetical probability density function underlying the multidimensional data is characterized by certain archetypal features. Recall that classical phenotyping is, in its essence, cluster analysis of data. Application of standard clustering methods assumes existence of lumps and peaks in the density function: therefore, clustering looks for a set of principal points  [ , ]. Dynamical phenotyping has a different basic assumption that the point density is characterized by the existence of continuous one-dimensional "ridges" which can diverge from or converge to each other in the data space.
However, the relation between 'biological' and 'chronological' age remains complex, especially in the pathological context [ ]. In our study we exploited the chronological patient age as any other clinical variable, and observed that indeed age correlates to some clinical trajectories but not to the others. Moreover, some clinical trajectories might be characterized by decreasing chronological age, which can be interpreted as aggravating clinical picture specific to younger patients. We can imaging other ways of using the age variable: for example, for learning the structure of the principal tree in a semi-supervised fashion. How to use the chronological age in the most informative way when analyzing both longitudinal and synchronic data remains an open question [ ]. Second approach to introduce physical time in the picture is using partially diachronic data as an additional annotation of a clinical dataset (the case of complete diachronic clinical data, representing longitudinal observations, is usually treated using a different and established set of approaches).
One source of information which can be relatively easily obtained is identifying pairs of data points corresponding to two subsequent states of the same patient, and recording the time lapse between two states. For example, part of the patients in a clinical dataset can be returning to hospital, with a previous record included in the dataset, so this information must be available. If the number of such pairs is sufficiently large then one can try to learn a monotonic function F k of pseudotime along each trajectory k, predicting the actual temporal label for each observation. Note that the connection between pseudotime and physical time can be different along different clinical trajectories. Moreover, the paired patient observation data can be used in the process of principal tree learning, by minimizing the number of paired patient observations belonging to distinct clinical trajectories.
Another limitation of the suggested approach is that the clinical trajectories are assumed to be diverging from some initial root state or states. In reality, convergence of clinical trajectories seems to be feasible (as in the case of the cellular trajectories). In this case, the model of the principal tree has to be generalized to some more general graph topologies (e. g. , existence of few loops). In case of ElPiGraph method, such modifications are easy to introduce technically: however, introducing graph structures more complex than trees requires careful consideration in order to avoid creating data approximators whose complexity will be comparable to the complexity of the data themselves [ ].
Potential implications Quantification of clinical trajectories represents the first step in using the concept of dynamic clinical phenotyping for diagnostics and prognosis. Predicting the probabilities of future clinical states for a particular patient together with their uncertainties, using the knowledge of clinical trajectories, can be a natural next step for future studies. These approaches can consider clinical trajectories as a coarse-grained reconstruction of the state transition graph for a dynamic system, described by, for example, continuous Markov chain equations. Some methodological ideas can be borrowed from the recent omics data studies [ ]. Recapitulating the multi-dimensional geometry of a clinical dataset in terms of clinical trajectories might open possibilities for more efficient applications of other methods, more oriented towards supervised machine learning.
For example, it can be potentially used for learning the optimal treatment policy, based on application of reinforcement learning as in [ ]. The existing large clinical datasets are frequently collected as a result of multi-site studies. In case of strong artefacts and biases caused by application of significantly different practices for data collection or other factors, specific methods of correction should be applied, integrated into the data analysis or even in the study design[ ]. However, dimensionality reduction methods based on averaging can in principle partially compensate for data heterogeneity if it can be modeled as a mixture of independent site effects that remain relatively small compared to the ranges of variable variations along the clinical trajectories. ElPiGraph in this respect has advantages over much more rigid Principal Component Analysis, being a nonlinear generalization of it for the case of datasets with complex geometries[ ].
Quantification of mixed type datasets Quantification of mixed type datasets, i. e. assigning a numerical value for nominal variables is a vast field where many solutions have been suggested [ ]. In the ClinTrajan package we used several popular ideas adapted to the aim of finding non-linear trajectories in the data. Firstly, for all non-binary categorical variables we suggest applying the "dummy" encoding (or "one-hot" encoding), e. g. introducing new columns, containing binary values one per each category (one of the categories might be dropped as redundant). Alternatively, if there is a sufficient number of numerical variables, Categorical Principal Component Analysis (CatPCA) can be applied  [ , ]. Secondly, for ordinal variables (including binary ones as particular case) we suggest to use either univariate or multivariate quantification.
For univariate quantification we assume that the ordinal values are obtained by binning a 'latent' numeric continuous variable possessing the standard normal distribution (zero mean and unit variance), following the approach described in [ ]. Let us consider an ordinal variable V which takes ordered values v < v. < v m , and each v i value has n i counts in the dataset. We quantify V by the values | x i = Φ - i- j= p j + p i , ( ) where p i = n i N , N is the total number of data points, and Φ(x) = √ π x -∞ e -x dx. If there exist many ordinal variables in the dataset, one can use multivariate methods for joint quantification of them. One of the most popular approaches is a particular variant of optimal scaling, aiming at maximizing the sum of squared pairwise correlations between all variables, including numerical and ordinal ones [ ].
ClinTrajan package includes its own implementation of this variant of optimal scaling which can be used to quantify ordinal variables in clinical datasets. The advantage of multivariate ordinal variable quantification with respect to the univariate one in that it can decrease the intrinsic dimensionality of the resulting data point cloud which can be beneficial for further applying of manifold learning methods, including the method of elastic principal graphs (ElPiGraph). The disadvantage of multivariate quantification of ordinal variables consists in necessity to have sufficiently large portion of data table rows without missing values. If this part is small then it might be not possible to quantify certain ordinal levels since they won't be represented in this complete part of the dataset, while this might still be possible with univariate quantification.
Imputing missing values in mixed type datasets The real-life clinical datasets are almost always only partially complete and contain missing values. Typically, these values are not distributed uniformly across rows and columns of the data matrix but rather form some non-random patterns, which can be even constructively used for the tasks of clinical data analysis [ ]. A typical pattern is existence of a column (or a row) containing abnormally large number of missing values. One can define two parameters δ row and δ column as maximally tolerable fraction of missing values in any row or column of the data matrix. The problem of finding the maximum size sub-matrix satisfying these constraints is not completely trivial but can be approximated by some simple iterative approaches. In practice, the trivial suboptimal solution consists in eliminating columns having the fraction of missing values larger than δ column , and then eliminating the rows having the fraction of missing values larger than δ row.
After constraining the maximum fraction of missing values in the data matrix, one can apply one of the available missing value imputation algorithms (imputers), which can be also classified into univariate and multivariate. For our purposes we advocate the use of multivariate imputers that allow us to avoid having strong data outliers destroying the manifold structure of the dataset. The standard scikit-learn collection provides two types of imputers: nearest neighbors imputation and iterative multivariate one, which can be in principle used for this purpose. In ClinTrajan package We add two alternative imputers based on application of Singular Value Decomposition (SVD) of order k. The first one which we will designate as "SVDComplete" is applicable if the number of rows in the data matrix having no missing values is sufficiently large (e.
The resulting elastic principal graphs are called elastic principal trees in this case. In ClinTranjan package we currently use only principal trees for quantifying trajectories and pseudotime, even though using more general graph topologies is possible. The advantages of limiting of the graph topology to trees are in that it is easy to layout the structure of the graph on a D plane and that any trajectory connecting two nodes of the graph is unique. The resulting explicit tree structure can be studied independently on the data. Also, an artbitrary vector x -not necessary belonging to the dataset X -can be projected onto the tree and receive a position in its intrinsic geodesic coordinates. The projection is achieved by finding the closest point on the principal graph as a piecewise linear manifold, composed of nodes and edges as linear segments connecting nodes.