TrendyGenes, a computational pipeline for the detection of literature trends in academia and drug discovery

Target identification and prioritisation are prominent first steps in modern drug discovery. Traditionally, individual scientists have used their expertise to manually interpret scientific literature and prioritise opportunities. However, increasing publication rates and the wider routine coverage of human genes by omic-scale research make it difficult to maintain meaningful overviews from which to identify promising new trends. Here we propose an automated yet flexible pipeline that identifies trends in the scientific corpus which align with the specific interests of a researcher and facilitate an initial prioritisation of opportunities. Using a procedure based on co-citation networks and machine learning, genes and diseases are first parsed from PubMed articles using a novel named entity recognition system together with publication date and supporting information. Then recurrent neural networks are trained to predict the publication dynamics of all human genes. For a user-defined therapeutic focus, genes generating more publications or citations are identified as high-interest targets. We also used topic detection routines to help understand why a gene is trendy and implement a system to propose the most prominent review articles for a potential target. This TrendyGenes pipeline detects emerging targets and pathways and provides a new way to explore the literature for individual researchers, pharmaceutical companies and funding agencies.


Gene annotation
We gathered the human gene synonyms from different sources (Ensembl, UniProt, HGCN, Entrez and Open-Targets; Fig. 1B) to sample the potential publications mentioning human gene names.Human genes had around 10 synonyms on average and many of those synonyms are ambiguous (Table 1): More than 30% of gene symbols had at least one promiscuous synonym, around 10% of the gene symbols are unsafe and have at least one gene synonym in the English dictionary, and almost 50% of gene symbols had a nested synonym.Combining these problems, almost 60% of the 19,082 gene symbols have one or more of these four types of ambiguity.To determine which synonyms are potentially ambiguous ("unsafe gene synonyms"; Fig. 1C) we did feature engineering to obtain variables that characterise unsafe synonyms (e.g.longer gene names are less probable to be ambiguous; (E) Prediction of per-gene publication trends using RNN.When a gene has significantly more publications or citations than expected by the model it is considered to be trendy.(F) Automatic topic detection of collections of publications.We used this algorithm to quantify the evolution of topics in trendy gene publications over time.(G) A review recommender system that uses information from the citation network and topic detection to recommend the most efficient set of reviews to explore the literature.
Table 2).Next, we used a positive-unlabelled bagging (PU) strategy following Mordelet et al. implementation 55 with a random forest classifier with the engineered features to calculate the probability of a gene synonym to be "unsafe" (see Methods).
To link every human gene to a subset of publications we implemented a disambiguation pipeline based on co-citation networks and machine learning (Fig. 1D).We gathered the titles, abstracts and keywords of the publications that had a match for any of the synonyms using regex with ElasticSearch (Fig. 1D).Nevertheless, this original set of publications potentially contains false positives: publications that contain an ambiguous gene synonym in their titles or abstracts, that do not refer to the gene of interest.
We assumed that true and false positives synonyms will tend to belong to different communities of publications from different research fields.To detect these communities we used co-citation networks (Fig. 1D): a weighted graph where the weight of the edges represents the frequency of two publications being cited simultaneously (co-cited) by a third publication.When two publications are repeatedly co-cited it strongly suggests that both belong to the same field of study 56 .We used the fast greedy modulation algorithm from iGraph to determine communities in the co-citation network and distinguished communities of publications focusing on the gene of interest by detecting the presence of "safe gene synonyms" in their titles and abstracts (Fig. 1D).The process is summarised in Fig. 2.
Finally, because we only used citations from open-access publications contained in PubMed Central (PMC) 57 , 46% of the publications were disconnected in the PubMed co-citation graph.To tackle this problem, we used again the inductive bagging positive-unlabelled approach to train multiple classifiers to associate the disconnected publications with the previously computed co-citation network components (Fig. 1C) using the words, phrases and one to four n-grams, contained in titles and abstracts.All available machine classifiers in Scikit Learn were used but logistic regression was selected due to its speed to accuracy ratio (Table 3).
To test the performance of the disambiguation pipeline we compared the disambiguation results with the gene-publication annotations from GeneRif 58 (manually curated annotations), DISEASES 59 (computational annotations), and UniProt 60 (computational and manually curated annotations) (Table 4).On average, the disambiguation recovers > 85% of all publications contained in these databases.Both GeneRif and Uniprot annotation do not necessarily contain a gene-synonym in the title or abstract, therefore those publications are out of our pipeline.Disambiguation results present on average a 70% precision with UniProt, the only collection of disambiguated publications of a similar magnitude.Finally, we included the disambiguated gene-publication annotations into the graph database.www.nature.com/scientificreports/

Trend detection
To detect incoming trends in the literature we gathered the publication dynamics of a given human gene from the disambiguated graph database (Fig. 1E).These time series include the number of publications, clinical trials, reviews and publications from big and medium-sized pharmaceutical companies, as well as, citations of publications coming from the mentioned categories per calendar year.Specifically, if a manuscript with author affiliations  to big pharma cites other publications these citations are categorized as big pharma citations.Conversely publications citing this manuscript whose authors are affiliated to big pharma are not categorized as big pharma citations.Time-series data from 1980 to 2013 was used to predict the per gene publication dynamics in each category between 2014 and 2019 using a Recurrent Neural Network model with an encoder-decoder architecture preceded by an attention layer, where both the encoder and decoder are composed of five hidden layers of Gated Recurrent Units (GRU).The time-series were created in a cumulative fashion, where each year contains the new publications and citations in addition to the previous ones.
For most genes, the model produces accurate predictions of the publication dynamics (Table 5), but for a small subset of genes the real number of publications or citations is significantly higher than expected (Fig. 3A).When the number of publications or citations exceeds the predictions, we interpret that the publication dynamics changed substantially in a way that cannot be explained simply by the gene's publication history, implying that a meaningful discovery in the field has recently occurred (Fig. 3A; orange).Trendiness is defined as the probability of the fold-change between predicted and real number of publications and citations for a given gene.We used this metric to identify the trendiest genes in the academic community-using all publications-, or in the pharmaceutical industry-using publications coming from pharmaceutical companies-(Table S1, supplementary material).
Finally, to identify trendy genes of pharmaceutical interest, we computed the normalised mutual information of genes and diseases in the titles and abstracts of publications (Fig. 3B).Disease names and their synonyms were obtained from the Medical Subject Headings (MeSH) ontology at the Bioportal 61 .MeSH ontology contains 4818 different disease nodes at different levels of the ontology.We created a dictionary for each disease with the preferred and alternative names (see Methods).The diseases were disambiguated in titles and abstract using the same disambiguation pipeline used with the genes.
We noticed that many trendy genes cluster forming trendy pathways when getting the gene-gene and genedisease association networks (Fig. 3C).We used enrichment of gene ontology (GO) terms for biological processes to uncover common pathways among the top 100 trendiest genes (Table S1, supplementary material).Among the most enriched GO terms in both academia and pharma are T cell co-stimulation, execution phase of necroptosis and pyroptosis.These biological processes are enriched in trendy genes which presumably reflect these fields of study are generating the most innovation and expectations in current biomedical research.

Topic detection
After the detection of gene trends, the next step was to understand why those genes might be trendy and curate possible mistakes in the disambiguation.With this aim we implemented a topic detection pipeline as an automatic, fast discovery tool to study groups of publications that mention the gene of interests (Fig. 1F).In this context, we used topic modelling algorithms.A topic is a collection of similar words, specific to a group www.nature.com/scientificreports/ of documents 62 .We used non-Negative Matrix Factorisation to generate a set of latent topics for each query (Fig. 4A; word clouds).We explored the evolution of the topics associated with some trendiest genes.For the immune checkpoint inhibitors (CD274, PDCD1, TGIT and CTLA4) the topic timeline suggests that there was a rapid decrease in the likelihood of publications discussing the biological role of these immune checkpoint inhibitors since 2010 (Fig. 4A in grey), which coincides with a notable increase in topics that discuss cancer therapies (Fig. 4A in orange) and monoclonal antibodies that target these four different transmembrane immunoglobulins (Fig. 4A in yellow).This way, the topic-detection pipeline is able to capture the evolution of the research from its biological description to the clinical application.
The topic timeline of the members of the necroptosis pathway (RIPK1, RIPK3 and MLKL; Fig. 4B) suggests that in the last decade there has been a decrease in the likelihood of publications discussing these genes in the context of apoptosis (Fig. 4B in grey), in favour of publications that discuss the newly discovered form of cell death, the necroptotic pathway (Fig. 4B in orange), as well as, the translational medicine perspective of this pathway as is suggested by words like mouse, treatment and activity or cancer (Fig. 4B, in blue).
Finally, the topic timeline the members of the pyroptosis pathway (CGAS, TMEM173, GSDMA and GSDMAD; Fig. 4C) shows a fast increase from 2013 of publications discussing the therapeutic opportunity in cancer immunotherapy with agonists for TMEM173 (Fig. 4B in grey), while again, the remaining topics seemed to contain information on the biochemistry and biological role of the genes.

Recommender system
In addition to the automatic topic detection, we designed a review recommender system to accelerate the screening of the publications that cover most of the information in a network (Fig. 1G).There are an average of 2.9 reviews citing any publication that mentions at least one gene name.The aim was to minimise the time reading and maximising the information within a gene subnetwork.The algorithm aggregates both topic and network information from the citation subgraph of the publications that mention the gene of interest to obtain the most query-centric reviews.The topic information comes from the latent topics obtained from the topic detection algorithm.The network information was captured by the PageRank scores of the subgraph (see Methods).This approach ensures that reviews citing publications with highest PageRank scores are prioritised.To further minimize the number of reviews for initial human analysis we avoid repetition of information by simultaneously maximising the cumulative PageRank score whilst minimising the overlap of combined citations.This way, we expect to obtain a small set of reviews that will cover the main topics and publications in the field.We used this recommender system to select the optimal subset of reviews to assess why genes might be trendy (see Discussion).An example of the output can be found for the genes in the discussion in the supplementary data file.

Discussion and conclusions
We present TrendyGenes as a first attempt to (i) establish a systematic analysis of contemporary topics associated to human genes and diseases, (ii) develop an alert system for emerging targets and trends in the scientific literature across the human, protein-coding genome, (iii) to use topic modelling to rapidly generate timelines of phrases that facilitate the understanding of why these genes are trendy.
We constructed a graph database containing PubMed data where publications are connected by citations and authors and are annotated with disambiguated human gene-names and diseases.We expect this new resource to provide new ways to navigate the scientific literature, detect and visualise networks of discussion and analyse networks of influence from key opinion leaders.Disambiguating author names from PubMed, MedRxiv, or BioRxiv would further improve the quality of the database.Machine or deep learning algorithms could be trained on already labelled data to improve on previously published approaches [63][64][65] and address this issue.
Similarly further improvements in gene-name disambiguation would assist precision and recall metrics on our validation set suffer for different reasons GeneRIF and DISEASES include fewer publications in comparison to the genome wide metrics identified in our pipeline and there will be a lot of potential "false positives".This makes the precision of our approach appear lower than what it may actually be.On the other hand, GeneRIF and Uniprot contain publications which either are not gene specific or do not mention the gene in the title or abstract.
However, the disambiguated genes and diseases can serve as labelled data for more sophisticated deep learning approaches to annotate biomedical entities.Gene and disease entities could be better annotated using both representation learning to capture the network topology and contextual information with transformer layers.Topic detection could be improved by using the state-of-the-art text summarisers with deep learning.
The number of publications per gene in aggregate is very predictable 66 .However, occasionally genes present significantly more publications than expected, meaning that a recent breakthrough occurred which cannot be accounted for from the publication dynamics.In this study, we show that trendiness can identify emerging targets from the literature for rapid profiling at genome-scale.We combined trendiness with gene-disease associations to prioritise potential drug targets: emergent genes associated with diseases but yet included in pharmaceutical publications are worthy of being investigated as potential targets.We observe that TrendyGenes usually cluster into the same biological pathways (Fig. 3C for CD274, PDCD1, CTLA4 and TIGIT).Here, using topic modelling and the recommendation system, we identify the trendiest genes and pathways and discuss some case studies to exemplify our pipeline.We selected genes pharmacological relevance by choosing genes with high trendiness both in the academia and the pharmaceutical industry with high association with disease and more than 100 publications.Reviews suggested by the recommender system for these genes are included in (whatReview2read.zip, supplementary material).CTLA4, PDCD1 (PD-1), CD274 (PD-L1) and TIGIT are among the trendiest genes in academia and pharma in 2019 (Fig. 3A).CTLA4, PDCD1, CD274 and TIGIT genes encode four different transmembrane immunoglobulins that act as co-inhibitory receptors: checkpoints or 'breaks' for the adaptive immune response that prevent T cells from exerting their functions 67,68 .CTLA4 competes with its analogous CD28 for CD80 and CD86 to prevent a premature activation of T cells 68 .PDCD1-CD274 interaction counters the positive signals that may have already activated T effector cells 68 .TIGIT interacts with CD155 to down-regulate natural killer cells and T lymphocytes 69 .Cancer cells attempt to impair these checkpoints and currently there are 7 FDA approved monoclonal antibodies that target three of proteins (CTLA4: Ipilimumab 70 ; PDCD1: Nivolumab 71 , Pembrolizumab 72 , Cemiplimab 73 ; CD274: Atezolizumab 74 , Avelumab 75 ) and multiple candidates targeting TIGIT (BGB-A1217 76 , OMP-313M32 77 , MTIG7192A 78 , AB154 79 ).Moreover, James Allison and Tasuku Honjo received the Nobel Prize in Medicine in 2018 for its research on immune checkpoint inhibitors 47 .

DNA sensing by cGAS-STING: cGAS, TMEM173, GSDMD, GSDMA
The cytosolic nucleic acid-sensing pathway leads to pyroptosis, a lytic pro-inflammatory type of cell death involved in antiviral, antibacterial, and anticancer response 107 .cGAS is a nucleotidyl-transferase that catalyzes production of cyclic GMP-AMP (cGAMP) upon the recognition of double-stranded DNA 107 .TMEM173 (STING) binds to cGAMP and promotes the activation of both TBK1 and IRF3, increasing the transcription of genes encoding type I interferons 107 .GSDMA and GSDMD are pore-forming effector proteins in the plasma membrane to release proinflammatory interleukins like IL-1β and IL-18 108 .The cGAS-STING pathway has been associated to multiple autoimmune and chronic inflammatory diseases like non-alcoholic fatty liver disease 109 , systemic lupus erythematosus 110 , vascular and pulmonary syndrome 111 , macular degeneration 112 , Bloom syndrome 113 , Aicardi-Goutières syndrome 114 , cancer 115 , DNA damage 116 , neurodegeneration 117 and beyond.Currently, there are ongoing clinical trials for TMEM173 [118][119][120] and GSDMD 121 although there are no reported trials for GSDMA nor cGAS.

Necroptosis: RIPK1, RIPK3, and MLKL
RIPK1, RIPK3 and MLKL form part of the tumour necrosis factor-induced necroptosis pathway [122][123][124] .This pathway has been associated with multiple pathologies: systemic inflammatory response syndrome 125,126 , ulcerative colitis 127,128 , psoriasis 128 , rheumatoid arthritis 128 , neurodegenerative diseases 129 and even cancer [130][131][132] .TNFR1, FasL, TRAIL, and TLR can all activate RIPK1 to decide the cell's fate: inflammation, apoptosis or necrosis 133 .If caspase-8 is inhibited, RIPK1 and RIPK3 form the necrosome that subsequently phosphorylates MLKL 134 .MLKL forms homo-trimers 135,136 , migrates to the plasma membrane 135,136 , binds to highly phosphorylated inositol phosphates 137 , creates pores in the membrane 138 and disrupts the cell integrity.The discovery of RIPK1 dates back to 1995 139 .Since then, four inhibitor programs have progressed through human phase II safety trials [140][141][142][143] .The first publication mentioning MLKL is more recent 144 and, despite the lack of kinase activity, pharmaceutical companies have cited its publications by 60 times more since 2013.Although there are no clinical trials yet, there are at least three known different chemical inhibitors 145 .

Mechanobiology: YAP1/WWTR1, PIEZO1 and PIEZO2
Cells use mechanical cues from their environment to guide behaviours such as proliferation and migration.Forces act as signals which are transduced to the nucleus where they control gene expression 146 .Mechanical forces are critical regulators of organ and tissue homeostasis, morphogenesis and regeneration, and are important aspects of diseases like cancer, metastasis, fibrosis and cardiac hypertrophy.YAP1/WWTR1 (TAZ) are transcriptional co-activators and mechanotransducers 147 .YAP/TAZ is hyperactivated in cancers 148 , its inhibition reduces atherogenesis 149 and fibrosis 150 , it triggers pulmonary hypertension 151 , and it is necessary for epithelial regeneration in the intestine 152 .PIEZO1 and PIEZO2 are two mechano-sensitive cation channels that play a key role in cell number regulation 153,154 and migration 155 , hearing 156 , neural 157 and vascular 158 development, somatosensory functions 159 , proprioception 160 and beyond.Piezo channels have been recently associated with multiple pathologies like arthrogryposis 161 , apnea 162 , congenital lymphatic dysplasia 163 , hyperalgesia 164,165 , malaria 166 , pancreatitis 167 , xerocytosis 168 , Gordon syndrome, Marden-Walker Syndrome, and Distal Arthrogryposis Type 5 169 .The discovery of mechanotransduction signalling pathways has received notable attention in the last years and may open the door to new therapeutic strategies to treat these diseases 147 .
Vol:.( 1234567890 www.nature.com/scientificreports/Trends in scientific literature are useful for pharmaceutical and biomedical companies.Moreover, this approach can offer crucial information to funding agencies to prioritise projects and a new way to study the research impact.Finally, individual researchers may benefit from a new methodology to explore the literature and from algorithms to maximise the efficiency of navigating over an increasingly vast biomedical literature.

Terminology
Here we use the term gene symbol to mean the approved symbol for any of the 19,084 human, protein-coding genes accepted by the HUGO Gene Nomenclature Committee.We refer to gene synonyms as any of the possible gene name variations by which the scientific community has ever referred to a given gene.Approved gene symbols are also included in the gene synonyms.For example: 'EGFR' is the approved gene symbol whereas 'EGFR' , 'Epidermal Growth Factor Receptor' , 'ERBB1' , 'ErbB-1' , 'c-erbB1' , 'HER1' , 'ERBB' are gene synonyms.We define promiscuous gene names as any gene synonym that is a synonym of more than one gene.This can include previous official gene symbols since these will not have been expunged from the literature.An example of this could be ' ARP1' which is a promiscuous gene synonym for the gene symbols 'NR2F2' , ' ACTR1A' , ' ACTR1B' , ' ANGPTL1' , ' APOBEC2' , ' ARFRP1' , 'PITX2' 47 .Unsafe gene synonyms are gene synonyms that may have a different meaning in other areas of research or in another context, for instance in standard English.The 'STAR' gene symbol is unsafe as opposed to its gene synonym 'Steroidogenic Acute Regulatory Protein' or CCP4 is both a gene synonym and the name for crystallography software.The final type of synonym we distinguish are Nested gene synonyms.These are gene synonyms that are part of another gene synonym.For instance 'insulin' is a nested gene synonym of 'insulin receptor' , 'TNF' is nested gene synonym of 'TNF Receptor Superfamily Member 1A' (gene symbol 'TNFRSF1A') and 'TNF Receptor Associated Factor 2' (gene symbol 'TRAF2').

Pubmed as a graph database
PubMed baseline 2020 53 comprises 30,419,056 publications for biomedical literature from MEDLINE and life science journals and 173,572,773 citations from the full-text archive of open-source publications PubMed Central (PMC).PubMed was imported into a graph database (Fig. 1A) for a fast performance in the retrieval of highly relational data like authorship and citation networks.In a graph database information is represented as nodes and edges, allowing the fast retrieval of queries about relationships.We loaded PubMed 2020 base-line into Neo4J 54 , an open source graph-database management system.We introduced four node types (publications, authors, human protein-coding genes, human diseases, medical subheadings), and four edge types (published, from authors to publications; cited by between publications, gene annotation from genes to publications; and disease annotation from diseases to publications).Furthermore, PUBLICATION nodes have the following attributes: PubMed identifier, title, abstract, affiliations, is_review, is_clinical_trial, big_pharma, med_pharma and date of publication.Profiling of the graph is included in Table 6, Database Profiling.Neo4J offers an interactive approach to navigate through PubMed (i) easily accessing references of publications, (ii) with the ability to query for specific genes and diseases already disambiguated, and (iii) with the aim of creating a knowledge graph for further exploration of gene-disease associations.The database is accessible to download at: https:// zenodo.org/ record/ 83626 79.

Gold standard sets
GeneRif 58 , UniProt, and DISEASES 59 were used as a golden-standard for validation.Table 6.Database profiling.We loaded PubMed 2019 base-line into Neo4J, an open source graph-database management system.We introduced four node types (PUBLICATION, AUTHOR, GENE, DISEASE), and four edge types (PUBLISHED, from AUTHOR to PUBLICATION; CITED_BY between PUBLICATION, GENE_PMID_ASSOCIATION from GENE to PUBLICATION; and DISEASE_PMID_ASSOCIATION from DISEASE to PUBLICATION).Furthermore, PUBLICATION nodes have the following attributes: PMID, TITLE, ABSTRACT, AFFILIATIONS, IS_REVIEW, IS_CLINICAL_TRIAL, BIG_PHARMA, MEDIUM_ PHARMA and DATE.The database is accessible at: https:// mega.nz/ file/ 4E8Qj CaQ# oqtm7 jof-lsG7y Sget8 uakh7 m26bD Lo1Hr Pu3mt dAV8.www.nature.com/scientificreports/

Pharmaceutical companies
A list of organisation names was generated from Cortellis 170 .Organisations with more than 100 patents in Cortellis were considered 'big pharma' and 'mid pharma' otherwise.

Gene annotation
Search for publications A ElasticSearch API search engine was used to retrieve PubMed IDs of publications containing a gene or disease synonym in their title, abstract or keywords (Fig. 1D; GET PMIDs and GET corpus).These PubMed IDs were later used to retrieve the publications' attributes from Neo4J using Cypher language through its python driver 171 (Fig. 1E; GET trends).Regular expressions were used avoid nested name ambiguity with lookarounds and fuzzy matching to account for case and punctuation and letter case variations (e.g.'ErbB-1' , 'erbB1' , 'ERBB1' , 'ErbB 1').

Unsafe synonym detection
We used 19,082 protein-coding human genes annotated by HUGO Gene Nomenclature Committee (HGNC).Gene synonyms were obtained from Ensembl, HUGO, Entrez, UniProt and Open Targets (Fig. 1C).Gene synonyms which were identical to disease names contained in the Medical Subject Headings (MeSH) database were eliminated.This mainly occurs when genes are named after diseases that are associated with e.g.'Li Fraumeni syndrome' as a gene synonym for gene TP53 172 or 'Marfan syndrome' in 'FBN1' 173 .Gene synonyms were classified into "safe" or "unsafe" categories using a modified version of the positiveunlabelled learning with bootstrap-aggregating as implemented by Mordelet et Vert (PU-learning; Fig. 1C) 55 .PU learning is a form of semi-supervised learning which iteratively finds positive examples within a-priori unlabeled data.To build a binary classifier able to distinguish the unlabelled class (U) into unsafe (P, positive) and safe (N, Negative) we engineered a series of features (Table 2, Unsafe features) such as the combined frequency of the characters in a gene synonym (example: 'ZNF' will be safer than 'EDA' because 'Z' and 'F' characters are less frequent in PubMed corpus than 'E' , 'D' , and ' A') or the probability of a gene synonym given that other gene synonym appeared in the text (the probability of 'STAR' given 'Steroidogenic Acute Regulatory Protein' is high but the probability of 'Steroidogenic Acute Regulatory Protein' given 'STAR' is low because 'STAR' is more ambiguous).
The PU learning was run for five iterations with a random forest classifier.The pure positive class (unsafe) was constructed combining gene synonyms present in the Enchant English dictionary 174 , gene synonyms with less than three characters, and promiscuous gene synonyms.In an active learning fashion, after each iteration, the top 1000 examples with the highest probability of being unsafe were manually relabelled if they were wrongly classified.For example, true positive unsafe synonyms like gene families (e.g.'G protein coupled receptor'), phenotypes (e.g.'Williams Beuren Syndrome') and other biological entities (e.g.'Cell surface antigen') were included in the true positive set for the next iteration.False positives like 'thymopoietin' or 'tubulin alpha-1C chain' were included into a new true negative class for the remaining iterations.
After the five iterations, a gene synonym was considered unsafe if: (i) it is included in the English dictionary from Enchant, (ii) it is a word with less than three characters, (iii) if the predicted score for the random forest classifier was higher than 0.5, and (iv) it is a promiscuous gene synonym .

Community detection
We produced weighted, undirected co-citation networks from unweighted, directed citation networks (Fig. 1D, GET communities).Subsequently, connected components were broken into communities using the fast greedy modulation algorithm implementation in iGraph 175 .

Gene annotation
Communities in co-citation networks represent different areas in the scientific literature.We used this feature to disambiguate large groups of publications (Fig. 1D, LABEL communities).We labelled all the publications in a community with the gene symbol of interest if the ratio of publications mentioning at least one safe synonym with respect of publications that only mention unsafe synonyms was higher than 0.1%.
Nevertheless, the co-citation graph is incomplete because PMC only contains citations of open-source publications.Because of that, 46% of the publications were disconnected from the co-citation graph.Disconnected publications that mention a safe-synonym were automatically linked to the gene symbol of interest.The rest of disconnected publications were linked to gene of interest using PU bagging strategy 55 with a binary logistic regression classifier based on the words in the text corpus (keywords, titles and abstracts) of communities already linked to the gene of interest and discarded communities.
Each corpus was pre-processed by (i) removal of non-alphanumeric characters, (ii) tokenization or split by whitespace, (iii) deletion of stop words from NLTK 176 , (iv) lower case conversion, (v) deletion of tokens whose length is less than three characters, (vi) deletion of token representing integers and (vii) stemming ('disambiguated' , 'disambiguations' , 'disambiguating' is converted to 'disambiguat').List of tokens (uni-, bi-, tri-, tetra-grams) with at least 2 counts and a frequency lower than 0.6 in the complete corpus were vectorised using TF IDF 177 .When there were less than 1000 unlabeled publications in the training set for the gene of interest, we generated an auxiliary negative class to augment the negative examples in the training data.This auxiliary negative class comprised a random sample of 1000 publications that mentioned genes different from the gene of interest.
These vectors were fed to all available machine learning classifiers from the Python library sklearn: Extra Tree Classifier, Gaussian Process Classifier, K-Nearest Neighbour, Logistic Regression, Ridge Classifier, Random Forest Classifier, and Support Vector Machine.All classifiers were trained with hyper-parameter tuning and threefold www.nature.com/scientificreports/cross-validation to avoid over-fitting in each of the 50 PU-bagging iterations (Table 7).The implementation of this PU learning algorithm is the same as the inductive bagging positive-unlabelled learning with bootstrapaggregating approach described by Mordelet 55 (PU-learning; Fig. 1D) and also the same as in the section of "Unsafe Synonym Detection" of the Methods (PU-learning; Fig. 1C).To maximise specificity and sensitivity simultaneously we selected the models with highest weighted F1 score (Scikit Learn) to maximise precision and recall at the same time.We selected the logistic regression (LOG, Table 3) classifier for the disambiguation pipeline given its accuracy-speed balance (Table 3).

Disease annotation
The same procedure used for gene entity recognition was used to detect disease entities, co-citation networks and machine learning.The Medical Subject Headings (MeSH) ontology was downloaded by querying their Rest-API available at BioOntology 178 .We constructed a key-value dictionary.Each disease was a node in the ontology.The disease synonyms were obtained from the 'Concept List Terms' field in the ontology to gather the preferred and the alternative ways of denoting the disease.We generated more synonyms of the diseases by reversing the order of synonyms with commas: 'Insipidus, Diabetes' to 'Diabetes Insipidus' .

Gene-gene-disease co-occurrence
Co-occurrence Co-occurrence of genes and diseases was computed using the simultaneous occurrence of gene/disease tags in publications after disambiguation, normalized by the total number of publications presenting those tags.We also computed mutual information metrics for gene-gene and gene-disease associations.

Gene mesh-parent association
Every disease MeSH term was associated with its lowest ancestor in the MeSH ontology under the node Disease 179 .After computing the gene-disease co-occurrence, each gene is linked with the most frequent ancestor disease term (Fig. 3B).

Trendiness
In this paper the trendiness for a gene is defined as the probability of observing the magnitude of fold-change between the predicted and the real number of publications for that given gene.The error in the predictions is inevitably higher with genes associated with small numbers of publications.To correct for this, we generated five bins based on the initial number of publications (percentiles 20, 40, 60, 80 and 100).We computed the distribution of the fold-changes between the predictions and observed reality in each of the five bins using a gaussian kernel density estimator available at Scikit Learn (bandwidth = 0.1, remaining parameters with default values).The area under the obtained probability density function is equal to 1. Trendiness is the area of the right tail of the probability density function bounded to the left by the observed fold change.This gives us an estimate of how extreme the fold change was for that gene in a specific bin.

RNN model
The model consists of an encoder-decoder architecture preceded by an attention layer.Both the encoder and decoder were composed of five hidden layers of Gated Recurrent Units (GRU).The model was implemented in Keras using the Tensorflow-GPU backend.Min-max normalization was used to rescale the time series before training.The optimizer was RMSprop and the loss was computed as the log error.30 percent of the time series was reserved for validation during the training.

Model optimisation
Input data was in both forms: cumulative and differential.In the cumulative model each year contains all the publication and citations published until then, while the differential model only contains the publications published in that year.Multiple normalisations were used ('none' , 'minmax' , 'log' , 'standard' , and combinations of them).Similar results were obtained with different normalisations and minmax was finally selected.Multiple Recurrent Neural Networks (RNNs) architectures were used (GRU, LSTM) in the form of encoder-decoder, with different numbers of neurons (1, 5, 10, 20, 50; Table 8).Models were compared with the Mean Accuracy Scaled Error (MASE), an unbiased method to compare time-series prediction models by comparing how much each model out-performs a naive model that repeats the last value.The 5-neuron-GRU with cumulative time-series was selected because it was the most parsimonious model with the smallest MASE.

Gene ontology terms enrichment
For the enrichment of gene ontology terms (Biological Process) associated with the 100 trendiest genes in academia (all publications) and the pharmaceutical industry we used the online tool GeneCodis 4.0 180 with default parameters.

Topic detection
We implemented algorithms to detect topics in collections of publications.This is useful to determine in which areas trendy genes are relevant.Furthermore, topic detection allows the fast identification of errors during the disambiguation.We used two different topic detection algorithms: Latent Dirichlet Allocation 62 (LDA), and

Figure 1 .
Figure 1.Workflow.Chart summarising the process from the downloading of the data to the detection and analysis of trends in the literature.(A) Creation of a graph database with the information contained in PubMed baseline 2020.(B) Acquisition of a comprehensive collection of human coding gene names and synonyms.(C) Automatic determination of potential ambiguous (unsafe) gene names.(D) Annotation of the graph database with unambiguous gene symbols by combining co-citation network topology and binary classifiers.(E)Prediction of per-gene publication trends using RNN.When a gene has significantly more publications or citations than expected by the model it is considered to be trendy.(F) Automatic topic detection of collections of publications.We used this algorithm to quantify the evolution of topics in trendy gene publications over time.(G) A review recommender system that uses information from the citation network and topic detection to recommend the most efficient set of reviews to explore the literature.

Figure 2 .
Figure 2. Disambiguation pipeline.(A) Citation network for a subset of PubMed IDs mentioning any of the gene synonyms of the gene symbol LRWD1, including ORCA.(B) Co-citation network of the same subset of PubMed IDs as in (A).(C) Communities for the co-citation graph obtained after using iGraphs fast greedy algorithm: killer whale community, orca plant cluster, LRWD1 in drosophila and LRWD1 in heterochromatin.(D) Number of safe synonyms per PubMed ID in title or abstract in the same co-citation network.(E) Citation network with reviews citing any of the PubMed IDs.(E) Review information as defined by the recommender system scaled from 0 to 1.

Figure 3 .
Figure 3. Trend detection and gene-gene-disease co-occurrence.(A) Logarithmic scatter plots showing the predicted number of publications, reviews, citations and citations from big pharma companies against real data in the year 2018.(B) Trendiness (log2(predicted/real)) for genes associated with groups of diseases (MeSH parent categories).Left; Average trendiness of publications, reviews, citations and citations from reviews.Right; Average trendiness of citations coming from big and medium sized pharmaceutical companies.(C) Gene-Gene-Disease co-occurrence network of the first neighbours of CD274.Orange nodes are diseases, grey nodes are genes and the size of gene nodes represents the trendiness The grey edges are gene-disease association, the blue edges are gene-diseases with the width of the edges reflecting the number of co-occurrences.

Figure 4 .
Figure 4. Topic time-lines.Topic time lines.Topic timelines for publications mentioning any of the genes for the immune checkpoint inhibitor (A), necroptosis (B) and pyroptosis (C) pathways.The x-axis represents the time in years and the y-axis represents the likelihood of a given topic.Colors represent different topics defined by the keywords contained in the correspondent word clouds.The latent four topics were obtained using Non-Negative-Factorization all publications annotated with the genes after disambiguation.Word clouds were created using the phrases with highest TFIDF for groups of publications belonging to each topic.All timelines show at least one rising topic after 2013 that represents the reason why these genes became trendy, their implications in human disease: immune checkpoint inhibitors and monoclonal antibodies (yellow and orange in A), activation of necroptosis (orange in B), agonists of STING1 in cancer (black in C).

Table 1 .
Gene synonyms are ambiguous.Manually discarded synonyms were labelled as unsafe during the unsafe gene synonym detection in an active learning fashion (see Methods).Unsafe aggregates the data from all the other categories.Data for 19,082 gene symbols and 185,549 gene synonyms.The total counts represent the number of individual synonyms when grouped by gene symbol and gene synonym.Promiscuous synonyms are counted as many times as they act a synonym.

Table 2 .
Unsafe features.Engineered features to evaluate the probability of a given gene symbol of being ambiguous (unsafe).

Table 3 .
Classifier comparison.Performance metrics for the 8 classifiers (Extra Trees Classifier, ETC; Gaussian

Table 4 .
Comparison of disambiguation methods.Average recall and precision of the disambiguation of our disambiguation with other databases.Low precision values for DISEASES and GeneRIFs are due to the smaller size of these databases.

Table 5 .
Performance of the predictions.Performance in the predictions of the publication dynamics.The model predicts the publications dynamics per gene between 2014 to 2019 using data from 1980 to 2013.Numbers represent the median 13,380 human genes.MASE mean accuracy scaled error, RMSE root mean square error, Total number of elements in the database up to 2013.Vol:.(1234567890)Scientific Reports | (2021) 11:15747 | https://doi.org/10.1038/s41598-021-94897-9

Table 7 .
Hyper-parameters for model selection.Hyper-parameters used during the model selection for Table3.ClassifierName column contains the name of the classifier and our acronym.Some classifiers were grouped since they have similar hyper-parameters like ExtraTreesClassifier (ETC) and Random Forest Classifier (RFC); and RidgeClassifier (RDC) and LogisticRegression (LOG).HParamName contains the names of the hyper-parameter names in the same format as sklearn version 0.24.2 (stable).ValuesToUse column contains the list of potential values of those hyper-parameters to be evaluated.Some values are specific for only one classifier and therefore have the acronym for the model in parenthesis (e.g.'lsqr' (RDC) and 'lbfgs' (LOG)).random_state and class_weight hyper-parameters were intended to have the same value across all models.The validation_fraction was used in MLPC to use the feature of early stopping : this created a sub-validation set under the training set different from the validation sets created for the cross-validation. hidden_layer_sizes[(10,),(50,),(100,),(10,10,),(50,50,),(100,100,),(10,10,10,),(50,50,50,),(100,100,100,)]

Table 8 .
Model size comparison.Global performance of the RNN model.Mean accuracy scaled error (MASE) is an unbiased method to compare time-series prediction models by comparing how much each model outperforms a naive model that repeats the last value.Values below 0 indicate that the method performs better than the naive model.Data for 13,380 human genes.