dataesr/scientific-paragraphs-categorization huggingface dataset:
The difference of the histological variants repartition in the BRAF V600E negative and positive tumor populations was assessed by a two-side Fisher's exact test which was very significant: pval < 0. 0001. Globally, similar results were obtained when comparing LNM positive N1 (n = 206) vs LNM negative N0 (n = 214) PTC samples and high risk (n = 24) vs low risk (n = 171) PTC samples (Additional file 1: figure  S2 ). As expected, the proportion of BRAF mutated tumors in the aggressive population is significantly higher than in the non-aggressive population (Fisher's exact test pval < 0. 01 in both comparisons). In order to compare the discrimination ability of related clinical and pathological information, we performed different PCA analyses on the TCGA dataset with the expression of the 14 down-regulated micro-RNAs.
Deep investigation of the miRNome revealed no major implication of isoforms variation or arm shift in PTC tumorigenesis and lymph node metastatic progression Since microRNA 5p-to-3p arm expression ratio modulations, NTA distribution and isomiRs distribution have been associated with tumorigenesis in several cancers  [23] [24] [25] [26] [27] [28] , we investigated the presence of such variations in our deepsequenced tumor and LNM samples (Table  1 ). We validated our results on independent samples from 14 patients by qRT-PCR when possible and with the TCGA dataset.
▪ 5p-to-3p arm ratio We calculated for each sample type and each pre-miR the expression ratio between both microRNAs (5p arm versus 3p arm) encoded by the same pre-miR. This calculation was only performed for microRNAs precursors encoding both a 5p and 3p mature microRNA. Expression of mature arms was performed in two ways: (1) using all distinct isoforms identified for a given microRNA (2) using only its canonical isoform reported in miRBase (v19). We found 9 pre-miR for which the 5p arm versus 3p arm expression ratio showed a fold change above 1. 5 or below 0. 66 between sample types (let-7d, mir-17, mir-27b, mir-92a-1, mir-181a-1, mir-181c, mir-221, mir-324, mir-1307), Fig. 7  Principal component analyses of whole TCGA dataset. a/b: the analysis was performed on the 59 normal samples and annotated tumors for BRAF mutational status (a) or histological subtypes (b).
Non annotated samples for one of these parameters were not used in the corresponding PCA. No additional filtering criterion was used. The analyses were performed with the expression of the 14 validated down-regulated microRNAs in PTC. Numbers in brackets correspond to the number of available samples in the TCGA dataset for each tissue type but only one of them (pre-miR-17) had an inversion of the dominant arm. We performed Taqman and SYBR Green (Exiqon) qRT-PCR validation experiments for 3 selected pre-miR (mir-324, mir-17, mir-1307) on the deepsequenced samples and the validation sample set. These two methods have been designed to target the canonical isoform reported in miRBase, unlike deep-sequencing. Although we confirmed the significant modulations observed for the 3 pre-miR in deep-sequenced samples, we only succeeded in validating the significant modulation observed for pre-miR-324 in the validation sample set.
However, neither pre-miR-324 nor the 8 other pre-miR had a significant modulation of the 5p-to-3p arm expression ratio in the TCGA validation dataset with or without filtering out samples with low tumor cells content (Additional file 1: Figure  S4 ). These data suggest that 5p-to-3p arm modulation is not a common feature of PTC tumorigenesis and lymph node metastatic progression.
b: only additions accounting for 1 % or more of all NTA in a given sample are represented. Other less frequent additions are grouped in the "other" category. Numbers correspond to the contribution, in percentage, of each type of addition averaged across all microRNAs having a minimum total mapped reads of 100. c/d/e: only isomiRs with a relative contribution to the total mapped reads of a mature microRNA of 1 % in at least one sample were considered. The number of isoforms per mature microRNA was then averaged by tissue type. The distributions were obtained in the 9 samples of the first sample set (c), in the TCGA samples dataset filtered based on clinical and pathological criteria (d) and in the complete TCGA dataset without additional filtering criteria (e). Numbers in brackets correspond to the number of available samples in the TCGA dataset for each tissue type.
We also screened NTA to identify specific microRNAs harboring different patterns of additions when comparing normal tissues to matching tumors and metastases, but we did not find variations of adenylation or uridylation accounting for at least 20 % of all mapped reads. This suggests that nontemplated additions in human are tissue and microRNA specific, as previously suggested  [37] [38] [39] [40] , but do not seem to be strongly related to tumor initiation or progression in the case of thyrocytes. Since the TCGA miRSeq public data only include per-sample isoform read counts and normalized expression, it was not possible to use them to reproduce this calculation (no information was available regarding the sequence of mapped reads).
First, by qRT-PCR on independent samples from 14 other patients. These samples were also carefully selected based on their normal or tumor cell composition. Secondly, we validated the results by an in silico analysis of the TCGA database which is the largest available dataset of PTC to our knowledge  [30]. qRT-PCR and in silico data were not always concordant. This can be due to different reasons. First, the number of samples analyzed is quite different between the 2 datasets. This variation changes the power of the associated statistics. Secondly, the methodologies are conceptually different  [21]. Finally, unlike qRT-PCR dataset, the TCGA dataset does not allow to adjust for interindividual variations by patient matched analyses without a great loss of information regarding LNM samples and not all the samples have been selected on the basis of their cell composition.
The implication of microRNAs in lymph node metastasis formation is only poorly characterized. It has been shown in different cancer types that microRNA expressions in LNM reflect the microRNA expressions in the associated primary tumors  [58] [59] [60]. Some studies suggest that LNM formation does not require major microRNA modulations  [60]  whereas another reports differentially expressed microRNAs between primary tumors and associated LNM  [58] , suggesting an implication of microRNAs in LNM formation. Concerning PTC metastasis development, very little is known. Our small RNA deep-sequencing results revealed similar microRNA expression profiles between primary tumors and associated LNM. In addition, we failed to validate by qRT-PCR 6 of the 7 selected microRNA modulations between tumors and LNM and only 2 (miR-873-5p, miR-876-5p) were validated in the TCGA dataset.
It is thus possible that a given isoform "A" reported in the TCGA data actually corresponds to isoform "B" with one or several non-templated additions in our own dataset. This limits the accuracy of comparisons of isomiR distribution. Nonetheless, we showed that the large majority of microRNAs had a level of adenylation or uridylation below 30 % of all mapped reads and it is rare to find microRNA with a high level of adenylation and uridylation at the same time. This is in accordance with previous results  [37, 66]. This particular issue only concerned isomiRs characterization and had no effect on expression profiles, obtained using all isoforms identified for a given mature microRNA. Additionally, qRT-PCR is inefficient for the validation of the results of the NTA and the isomiR analyses.
However, to our knowledge, this is the first work which explores 5p-to-3p expression ratios and NTA distribution in PTC. IsomiR distribution has already been studied by Swierniak et al. who showed that this distribution changes during PTC tumorigenesis  [27]. In accordance with them, we showed that the vast majority of mature microRNAs, covered by a minimum of 100 reads, presented between 2 and 8 isomiRs. However as opposed to what they report, we did not find modulation of the distribution between sample types, even when considering all the samples from the TCGA dataset. This discordance could be explained by the number of samples considered and by technical differences. Indeed, Swierniak et al used normal samples, normal tissues adjacent to tumor samples and tumors from 14 patients. In addition, they required perfect matches without any 3' end NTA trimming when aligning reads.
Acknowledgements Authors want to thank  Professor J. E Dumont  for his conceptual support,  Bernadette Bournonville  and  Chantal Degraef  for her excellent technical support. This work was supported by the  "Fonds de la Recherche Scientifique Médicale" (FRSM)  ;  Fonds de la Recherche Scientifique-FNRS  (grant number  19556102 ) ;  WELBIO  ;  Plan Cancer Belgique  ;  Télévie ; Fonds pour la formation à la Recherche dans l'Industrie et dans l'Agriculture  ;  Fondation Rose et Jean Hoguet and les Amis de l'Institut Bordet.
Introduction The hypothesis of modularity in muscle activity relies on the premise that the central nervous system (CNS) stores a limited set of modules and recruits them to execute the motor tasks at hand  (Bizzi et al. , 2008; Tresch and Jarc, 2009; Bizzi and Cheung, 2013). Putative modules are typically identified by means of dimensionality reduction algorithms applied to electromyographic  (EMG)  recordings  (Tresch et al. , 2006). This procedure yields modular decompositions that approximate the EMG data as accurately as possible despite the significant measurement noise that affects such data. Different types of modularity have been proposed assuming the existence of either spatial  (Tresch et al. , 1999; Ting and Macpherson, 2005) , temporal  (Ivanenko et al. , 2004 (Ivanenko et al.
In total, the experimental protocol specified four targets × two directions, that is, eight distinct movement directions. Two speeds were tested, so that the 16 distinct motor tasks were denoted by n 1 ,. , n 8 and f 1 ,. , f 8 for normal and fast speeds, respectively. Each motor task was composed of 40 trials. Thus, we had a number of muscles M = 9 and a number of samples S = 16×40 = 640, for each participant. Here, we illustrate our method on the dataset recorded from one representative subject and then test it on all other subjects to assess its robustness. FIGURE 1 | Schematic illustration of the steps of the proposed module identification approach. Left: The activity of several muscles is recorded simultaneously during performance of several repetitions of a variety of motor tasks. For illustration, here we show the activation of nine muscles when performing three tasks (task1 in red, task2 in blue, task3 in green).
To discriminate between tasks, the dispersion of module activations within each task (S (1) w , S (2) w ,S (3) w shown as red, blue, and green ellipses, respectively) should be minimized, whereas the dispersion of module activations across tasks (S b shown as a dashed-line black ellipse) should be maximized. E 2 LDA is the task discrimination cost expressed as a linear combination of S w and S b. Right: The proposed method integrates the data approximation and the task discrimination objective to identify invariant spatial/temporal modules and activation coefficients that combine them in single samples. For illustration, we show two spatial and two temporal modules combined with 2 × 2 = 4 activation coefficients to reconstruct each sample of muscle activity. Body kinematics was recorded by means of a Vicon (Oxford, UK) motion capture system.
Six passive markers were placed on the finger tip, wrist (over the styloid process of the ulna), elbow (over the lateral epicondyle), right shoulder (on the lateral epicondyle of the humerus), back of the neck, and left shoulder. The kinematics data were low-pass filtered (Butterworth filter, cut-off frequency of 20 Hz) and numerically differentiated to compute tangential velocity and acceleration of the fingertip. Movement onset and movement end were identified as the times in which the fingertip velocity profile superseded 5% of its maximum. The mean movement duration varied across subjects from 370 to 560 ms. Electromyographic activity was recorded by means of an Aurion (Milan, Italy) ZeroWire wireless surface electromyographic system. The EMG signals were recorded from the following muscles: 1.
finger extensors of the posterior forearm (extensor digitorum) (Fe), 2. brachioradialis (Br), 3. biceps brachii (Bi), 4. triceps medial (Tm), 5. triceps lateral (Tl), 6. anterior deltoid (Ad), 7. posterior deltoid (Pd), 8. pectoralis (Pe), 9. latissimus dorsi (Ld). EMG signals were digitized, amplified (20-Hz high-pass and 450-Hz low-pass filters), and sampled at 1000 Hz (synchronized with kinematic sampling). Subsequently, in order to extract the signal envelopes, the EMGs for each sample were digitally full-wave rectified, low-pass filtered (Butterworth filter, cut-off frequency of 3 Hz, zero-phase distortion;  Ivanenko et al. , 2004)  and their duration was normalized to 50 time steps. The data were then normalized in amplitude on a muscle-permuscle basis by dividing each single-sample muscle signal by its maximal value attained throughout the experiment.
NMF-based module extraction To uncover the spatial and temporal modules, we used an NMF-based algorithm, called sNM3F  (Delis et al. , 2014). The sNM3F algorithm takes as input the preprocessed EMG matrix and the desired number of modules (N and P) and gives as output the matrices W, A s , W representing the spatial modules, sample-dependent activation coefficients, and temporal modules, respectively. The algorithm aims to reconstruct the recorded EMG data by minimizing the following NMF-based error, while imposing non-negativity constraints on W, A s , W: E 2 NMF = S s=1 ||M s -WA s W|| 2 fro. (3) where ||. || fro is the Frobenius norm. The metric typically used to assess how accurately the modular decomposition reconstructs the original data is the Variance Accounted For (VAF) (d  'Avella et al.
(2013a Delis et al. ( , 2014)). The core idea of the method is that decoding performance should significantly increase only if inclusion of an additional module describes reliably some task-related EMG variations not described by other already included modules. Hence, this procedure ensured the detection of modules that explain only the "task-relevant" variability and the exclusion of other sources of noise that produce "task-irrelevant" variability. This twostep process (NMF-based plus LDA-based) has been shown to effectively identify spatial and temporal modules that explain task differences in muscle activations  (Delis et al. , 2014). For the dataset under investigation, P = 3 temporal modules and N = 4 spatial modules were shown to be the optimal choice (see Results, Spatial and Temporal Modules).
Alternating constrained least-squares (ALS) We consider an alternative formulation to specifically preserve the normalization constraint during the update of W and W. To this aim, we considered a constrained optimization problem. min W ||M -WR|| 2 fro s. t. T i=1 Wij = 1 ∀ j ∈ [1, P] (10) After vectorization using Kronecker products, this can be rewritten as: 2RR ⊤ ⊗ I T I P ⊗ 1 T I P ⊗ 1 ⊤ T 0 vec( W) λ = vec(2MR ⊤ ) 1 P (11) where 1 P is a vector of P ones and ⊗ denotes the Kronecker product. This is a standard system of linear equations for which a leastsquare solution can be efficiently obtained (e. g. , mldivide Matlab function). Similarly, to update W, one can solve the following leastsquare problem: I N ⊗ 2Q ′⊤ Q ′ 1 M ⊗ I N 1 ⊤ M ⊗ I N 0 vec(W) λ = vec(2Q ′⊤ M ′ ) 1 N. (12) Since non-negativity is not ensured after such a process, it can be imposed a posteriori by setting all negative elements to zero.
Non-Linear Programming (NLP) To derive our third method, we formulated the problem as a single large-scale constrained optimization problem. We implemented the gradients of the cost function and of the non-linear constraints to perform the actual optimization using fmincon Matlab function with calculation of the Hessian by a limited-memory, large-scale quasi-Newton approximation (namely, the interior-point algorithm with lbfgs method). The NLP problem was thus formulated as: min z E 2 (z) s. t. z i ≥ 0 ∀ i ∈ [1, TP 2 N 2 SM] and φ(z) = 0. ( 14 ) where z T vec( W) ⊤ , vec(A) ⊤ , vec(W) ⊤ , i. e. , the concatenation of the vectorizations of W, all the A s , and W. The function φ encoded the additional constraints to normalize the rows of W and columns of W so that they summed to one (linear constraint).
Comparisons of the three methods revealed that they were effective and yielded the same modular decompositions. The best choice of algorithm implementation (in terms of speed and robustness) may depend on the data set, but overall the MULT algorithm yielded good performance. Matlab code implementing DsNM3F including these three optimization algorithms is available at the last author webpage (http://hebergement. u-psud. fr/berret/software/ DsNM3F. zip).
Comparison between Decompositions We defined VDM = VAF × DEC (0 ≤ VDM ≤ 1), a measure of the goodness of an EMG decomposition that incorporates the VAF (data approximation) and the decoding performance (task discrimination). Higher values of VDM correspond to more plausible decompositions that achieve the best trade-off between data approximation and task discrimination. Hence, we used this metric to select the best output of DsNM3F with respect to parameters γ and δ, and compare it to standard space-by-time decomposition identified by sNM3F. To assess the similarity between modules identified by DsNM3F and sNM3F, we used as similarity index the correlation coefficient between pairs of modules (either spatial or temporal). To compare activation coefficients, we computed the average coefficients for each of the 16 tasks, which gave a 16-dimensional vector for each decomposition and computed the correlation coefficients for pairs of such vectors.
To assess the power of a modular decomposition in terms of task prediction, we used Support Vector Machine (SVM) classification (implemented in Matlab using the LIBSVM toolbox  Chang and Lin, 2011). We applied linear SVM to the activation coefficients combined with a repeated random subsampling cross-validation procedure  (Cortes and Vapnik, 1995; Theodoridis and Koutroumbas, 2008). More precisely, we split the dataset into training and test data 20 times making sure that both training and test data had always the same task proportions. More precisely, we tested the cases of 5, 10, and 30 training samples per task. For each one of the 20 splits, first we applied the module identification algorithms (sNM3F or DsNM3F) to the training data. Then, we kept only the extracted modules since they represented the basic modules on which the muscle patterns could be projected.
Therefore, we projected all the data (both training and test) on them without taking into account task information (as it was anyway unavailable for the test data). To do this, we used the sNM3F with fixed spatial and temporal modules, i. e. , updating only the activation coefficients. This procedure gave us new activation coefficients for both the training and test data sets. We used the activation coefficients of the training dataset to train a linear SVM classifier and the ones of the test dataset to compute the classification performance (percentage of correctly classified samples). We repeated this process 20 times and computed the average (±SD) classification score (denoted by CLS).
Results To illustrate the newly developed DsNM3F algorithm and have a better understanding of its functionality, first we tested it on simulated data for which the modules and their activations are known by construction. Then, we applied to the complete EMG dataset recorded during performance of 16 arm pointing tasks (40 trials per task), which has been analyzed before using the space-by-time decomposition model  (Delis et al. , 2014). Both datasets allowed us to make direct comparisons between the outputs of DsNM3F and sNM3F.
Trade-off between task discrimination and data approximation To exemplify the usefulness of the new algorithm, we generated an artificial EMG dataset consisting of two spatial modules and two temporal modules (Figures  2A, B ) and simulated their activations during performance of four motor tasks, each repeated ten times with relatively small inter-trial variability (Figure  2C , black bars) that allows perfect task discrimination. For simplicity, we modeled the spatial modules to consist of two muscles and the temporal modules to consist of two time points. The EMG signals (constructed by combining spatial and temporal modules using the single-trial activations) were corrupted with large additive and multiplicative noise (40% of the signal) to model different sources of measurement noise in the EMG recordings.
We applied the standard sNM3F and the new DsNM3F to the resulting EMG datasets and found that both algorithms reconstructed accurately the original spatial and temporal modules despite the presence of significant noise (Figure  2 , dark gray for DsNM3F, light gray for sNM3F, and black for original modules). DsNM3F approximated the original activation coefficients more accurately than sNM3F indicating that the insertion of the task discrimination objective resulted in a better reconstruction of the task-dependent structure of the module activations. This result is also reflected in the VAF and DEC values of the two decompositions. In light of the data generation, an imperfect VAF is expected even if both the modules and the true activation coefficients are correctly recovered because of the addition of unstructured noise to the EMG signals.
Influence of parameters on VAF and task decoding DsNM3F has two free parameters γ and δ that need to be carefully selected in order to achieve the optimal solution that balances task discrimination with data approximation. In our framework, we measure task discrimination as percent correct decoding performance (DEC) and data approximation as the VAF by the decomposition. As our goal is to achieve the highest possible task discrimination without a significant data approximation loss, we chose to maximize VDM, a metric that incorporates both objectives (see Materials and Methods for details). To determine the optimal pair of parameters, we developed an automated optimization procedure. First, we used a coarse grid in (δ, γ )-space to find an approximate value maximizing VDM. Second, we ran an optimization starting from the latter initial guess in order to tune more finely the parameters and try to get the best possible VDM value.
Interestingly, decoding shows a higher dependence on γ than δ indicating that decreasing the withintask variability plays a bigger role in decoding performance than increasing the between-task variability. The maximal VDM appears as a compromise between the two metrics and corresponds to γ = 0. 01 and δ = 0. 01 (illustrated with a star on the graph). These parameter values were actually very close to the global optimum (the best solution we found was for (γ , δ) = (0. 01, 0. 006) giving a VDM value very close to the one found using the coarse grid procedure), which finalizes our process for finely tuning the two free parameters. This approach was successful in determining the optimal values of γ and δ also for the data recorded from the remaining five subjects. Parameter values were found to be quite robust across datasets: (γ , δ) = (0.