Companion canines: an under-utilised model to aid in translating anti-metastatics to the clinic

Alongside shared environmental exposures and similar microbiomes [2], the spontaneity and heterogeneity of companion canine cancers arguably has more relevance to human patients than traditional laboratory models, such as rodents. Indeed, the ‘One Health Initiative’ recognises the link between the health of humans, animals, and the environment (http://www.onehe althi nitia tive.com/). The dog genome is well annotated and the human genome shares more ancestral sequence with the dog than the mouse [3], a reflection of which is that for many cancer-associated proteins there is a greater degree of sequence similarity between the human and canine homologs than between the human and mouse proteins [4]. At the protein level, antibodies used in humans for immunohistochemistry and flow cytometry often also work in dogs [5, 6], whereas mice frequently require an antibody specific for the mouse protein, further emphasizing the similarity between human and dog genes. Importantly, many canine metastatic tumours show striking aetiological, epidemiological, genetic and histological similarity to their human counterparts [7]. In addition, as is the case for human cancers, the presence of specific immune cell types in the tumour microenvironment of dogs has been linked with prognosis [8–10].

The dog as a model organism for translational human cancer research
Alongside shared environmental exposures and similar microbiomes [2], the spontaneity and heterogeneity of companion canine cancers arguably has more relevance to human patients than traditional laboratory models, such as rodents. Indeed, the 'One Health Initiative' recognises the link between the health of humans, animals, and the environment (http://www.onehe althi nitia tive.com/). The dog genome is well annotated and the human genome shares more ancestral sequence with the dog than the mouse [3], a reflection of which is that for many cancer-associated proteins there is a greater degree of sequence similarity between the human and canine homologs than between the human and mouse proteins [4]. At the protein level, antibodies used in humans for immunohistochemistry and flow cytometry often also work in dogs [5,6], whereas mice frequently require an antibody specific for the mouse protein, further emphasizing the similarity between human and dog genes. Importantly, many canine metastatic tumours show striking aetiological, epidemiological, genetic and histological similarity to their human counterparts [7]. In addition, as is the case for human cancers, the presence of specific immune cell types in the tumour microenvironment of dogs has been linked with prognosis [8][9][10].

Advantages of involving companion canines in the development of anti-metastatics
An important challenge in the development of anti-metastatic drugs is target identification. In order for a drug to be developed and trialled, there has to be a molecule worth targeting, and the definition of what constitutes a 'suitable' metastatic target is not widely agreed. However, comparison of the molecular differences between the metastatic tumour cell and the primary mass from which it evolved may elucidate potential therapeutic targets. In humans, it is often not possible to easily obtain both primary tumour and metastatic tumour from the same individual, and even more difficult to obtain the latter before any treatment has commenced. Tissue samples from the metastatic tumour as well as the primary tumour of the same animal are often available at postmortem, although whether or not they are collected depends on the institution and the resources they have available to sustain this as a routine practice. Humane euthanasia of dogs with terminal cancer and declining health is common, and when this takes place at a suitably resourced veterinary college with a pathology department many owners consent to necropsy (autopsy), realising that understanding the cause of death of their beloved pet can help with research trying to prevent the death of other pets.
One advantage of using companion canines in clinical trials of anti-metastatics relates to detecting the presence of metastases. The radiographic modalities for imaging metastases in humans (such as MRI, CT, PET-scan) are the same as that used in dogs and are commonly available at veterinary referral hospitals; these have a dual purpose, in being used to monitor metastatic burden as well as measure response to therapy. As an alternative to imaging-mediated treatment response monitoring, it has recently been demonstrated that dogs with osteosarcoma (OSA) that were free of metastases at diagnosis after having received limb amputation and adjuvant chemotherapy, can develop detectable levels of circulating tumour cells (flow cytometric detection of intracellular collagen 1) [11] in their blood before the development of overt metastases or death. The predictive benefits this offers means the relative success of a drug being trialled could be determined without having to wait for the metastatic lesion to become detectable by imaging and/or the dog to suffer clinical signs of tumour burden. In addition, identification of prognostic markers for metastasis development in dogs may also be translatable, facilitating monitoring of high risk human patients to enable early detection of metastasis.
The demonstration that gene expression signatures of metastatic potential could be detected in primary human tumours [12] lead to the identification of metastasis-associated gene expression signatures by comparison of gene expression in human primary tumours that did and did not metastasise [13,14]. Veterinary researchers have performed similar comparisons of metastasising and non-metastasising primary tumours to identify molecular genetic and epigenetic events associated with canine tumour metastasis, both as targets for anti-metastatic therapeutics and potential prognostic biomarkers [15][16][17][18][19]. Such studies will expedite the identification of potential targets for drugs aimed at preventing human tumour metastasis by strengthening the 'candidature', as a target, of genes that are associated with the metastasis of a given tumour in both humans and dogs.
As outlined above, matched primary tumour and terminal metastases from the same dog can sometimes be obtained. This allows for detailed examination of the metastatic tumour cells and the tumour microenvironment including the immune infiltration present, to gain some understanding why a therapy failed in some individuals but greatly prolonged disease-free survival in others. In pursuing the goal of treating metastatic lesions, we recognize that metastases are different from primary tumours, so we need to know the molecular details of metastases that resisted therapy, and this requires sampling these lesions. It is important to consider that dogs develop metastatic disease spontaneously and all the processes in the invasion-metastasis cascade occur naturally, and within the context of a competent immune system. This contrasts with most rodent models where intravenous, intracardiac, or intraosseous injections of already highly malignant and metastasis-competent tumour cells are used to model metastases. It stands to reason that studying spontaneously occurring metastases in dogs increases the likely relevance of preclinical testing of anti-metastatic agents towards treating metastatic progression in human patients. For human and canine metastatic cancers that have a shared metastasis-associated molecular drug target, and which metastasise to similar sites, the rationale for the preclinical testing of an anti-metastatic agent on canine patients is highly compelling.

Limitations of dogs in aiding cancer drug development
Whilst there are many reasons why dogs represent a good animal model of cancer in humans, there are a number of issues which continue to preclude the involvement of dogs in cancer drug development.
The identification of molecular drug targets that are shared between human and canine tumours will require the analysis of many more canine tumours. The number of canine tumours that have been subjected to genome-wide molecular profiling are but a small fraction of the number of human tumours that have been analysed.
Further along the drug development pathway, there has been some variability in the rigour of canine oncology clinical trials, and the outcome of such trials may not necessarily inform the outcome of a human clinical trial featuring the same tumour and therapeutic. Although much of the infrastructure of canine clinical trials, such as ethical approval, informed consent, data/sample management, and statistical analysis, is comparable to that of human clinical trials, historically there have been no standardized guidelines for conducting a canine clinical trial and thus no specific quality control and/or quality assurance procedures are required to be incorporated into all veterinary clinical studies. For example, trials with non-randomized cohorts are not uncommon, and may lead to skewed results. Non-randomized (e.g. phase II) clinical trials also occur in human medicine, however successful phase II human trials more often lead to randomized phase III trials compared to veterinary medicine, where historically such large trials are less common. As a result, there is an increasing awareness of the need to establish rigorous, standardized veterinary clinical trial designs to ensure valid results are obtained, and to this end workshops have been set up to establish 'Best Practice Recommendations' [20]. One example of this increased standardization and rigour is the Comparative Oncology Trials Consortium, which is part of The National Cancer Institute in the USA [21]. In addition, the American Veterinary Medical Association (AVMA) maintains a clinical trials website (https :// ebusi ness.avma.org/aahsd /study _searc h.aspx), which is similar to www.clini caltr ials.gov. It should be noted that both the AVMA and clinicaltrials.gov are essentially repositories that do not necessarily vouch for the rigour of the trials that are listed, and there is no 'standard' for inclusion of a trial in these lists. However, it is anticipated that this growing formalization of infrastructure for companion animal clinical trials will greatly increase the reliability of such trials and promote the widespread acceptance of companion canines as valuable models for oncology clinical trials with translatability to humans.
Whilst the study of a spontaneous, clinicopathologically-similar tumour in a dog is intuitively more relevant to a human tumour than the study of an induced 'model tumour' in a mouse, differences still exist between the species that can impact on the translation of clinical trial results from dogs to humans. For example, the acyclic nucleotide analog (GS-9219), which showed marked cytotoxic effects in human lymphoma and leukaemia cell lines in vitro, was evaluated in dogs with naturally occurring, advanced-stage lymphoma and found to display significant efficacy with an acceptable adverse-event profile [22]. However, in subsequent stage I and II clinical trials in humans with haematological malignancies the drug showed marked toxicity and development was stopped. Consequently, we need to accept that, as with any animal model, there will be occasions when differences between humans and dogs mean that drugs cannot always be shared.
The question is, acknowledging these limitations and differences that exist between the two species, do anti-metastatic clinical drug trials in companion canines represent a viable option that should be more vigorously pursued? Perhaps this question is best answered by looking at what we have learned from anti-metastatics clinical trials in companion canines to date.

Canine osteosarcoma: a highly relevant model for testing anti-metastatics
The majority (> 95%) of canines with appendicular osteosarcoma (OSA) have microscopic pulmonary metastasis at the time of diagnosis and a median overall survival of 4-6 months without adjuvant chemotherapy. Similar to humans, the treatment for OSA in dogs typically involves resection of the primary tumour (generally amputation) with chemotherapy, which increases median survival to 10-12 months. Unlike humans, OSA is a much more common tumour in dogs with an estimated > 10-fold higher incidence than in humans [23]. Fifty percent (14/28) of trials currently listed on the AVMA clinical trials website pertain to OSA, and it is these highly metastatic canine OSAs that have afforded the opportunity to rapidly progress the development of certain targeted therapeutics in human trials. For example, at the end of 2018, Advaxis licensed HER2-targeted agent ADXS31-164 (ADXS-HER2) to OS Therapies for development in human clinical trials for the treatment of patients with recurrent, completely resected OSA https ://advax is.com/clini cal-trial s-3/her2-assoc iated -cance rs/.
Toceranib phosphate (Palladia ™ ) is a tyrosine kinase inhibitor (TKI) that was approved by the FDA in 2009 as the first canine-specific cancer drug for treatment of cutaneous mast cell tumours. 1 [24] In a clinical trial at the Flint Animal Cancer Centre at Colorado State University, dogs with lung metastatic OSA treated with toceranib and high dose losartan (to suppress the activity of inflammatory monocytes 2 ) demonstrated a biological response rate of 50% and a measurable response rate of 30%, with acceptable toxicity (AVMA Animal Health Studies Database study: AAHSD000259). With this success, a second multi-institutional clinical trial was initiated in late 2018 using high-dose losartan and toceranib (AAHSD004794), and at the same time, in conjunction with paediatric oncologists at Children's Hospital Colorado, a clinical trial was being designed for the use of high-dose losartan and the TKI inhibitor, sunitinib, for paediatric bone cancer patients https ://www.csuan imalc ancer cente r.org/ blog/new-hope-for-canin es-and-kids-with-bone-cance r.
The apoptosis-promoting drug Procaspase-activating compound 1 (PAC-1) is another example of where success after rigorous evaluation in canine cancer patients paved the way for evaluation of the drug in human clinical trials. As a single agent, PAC-1 has shown considerable activity in canines with lymphoma. However, it also may potentiate other therapies. A PAC-1/doxorubicin combination treatment lead to a biologic response in 3/6 dogs with metastatic OSA, and 4/4 dogs with lymphoma [26], whilst a PAC-1/ temozolomide (TMZ) combination achieved biological responses in 3/3 dogs with glioma [27]. Although there are trials at different Institutes within the USA that are currently recruiting canines with lung metastatic OSA to further assess the effectiveness of the PAC-1/doxorubicin combination, Phase I clinical trials using the PAC-1/TMZ combination in humans with high grade gliomas (glioblastoma multiforme or anaplastic astrocytoma after progression following standard first line therapy) have already begun (ClinicalTrials. gov Identifier: NCT03332355). The hope is now that with such a successful track record, it may be possible to address the effectiveness of these combinations in clinical trials of patients with other types of metastatic cancers.

Current anti-metastatic clinical trials in dogs
Orally administered rapamycin (Sirolimus™), an mTOR inhibitor, is currently being investigated by the Comparative Oncology Trials Consortium, in many different Universities and centres across Canada and the USA, for dogs with amputation-confirmed OSA and no evidence of metastatic disease. The dogs receive standard chemotherapy (carboplatin chemotherapy), followed (or not) by oral rapamycin, with the primary outcome event being the time until development of metastasis (the trial has completed enrolment and the results are due to be released). The National Cancer Institute is also currently running a multi-centre trial (in Canada and the USA), to evaluate a recombinant, attenuated Listeria monocytogenes expressing a chimeric human HER2/neu protein as an adjuvant treatment for dogs with OSA, specifically assessing the anti-metastatic effects of the vaccine compared to standard treatment alone (Trial ID: COTC026).

Mutual benefits
Of course, the suggestion that companion canines should be used in clinical trials of anti-metastatics must not solely be for the benefit of humans-pet dogs are not just alternative mouse models! Dogs with cancer should be involved in the preclinical testing of an anti-metastatic agent where it has been demonstrated that the canine cancer shares the metastasis-associated drug target with the human cancer that the therapy is ultimately intended to treat. There is the possibility that whilst a drug may show promising results in companion canines, it may later be found to fail in human clinical trials, although one would predict that the overall failure rate of this approach in humans would be less than what is currently seen with conventional pre-clinical models. However, most pet owners would probably agree that developing a new treatment that was ultimately found to only benefit dogs would still be a valuable outcome.

What is likely to promote increased use of companion canines in anti-metastatics drug development?
A key factor to facilitate increased use of companion canines in anti-metastatic drug development is funding. Funding veterinary research for companion animals can be challenging. It is important to appreciate that whilst they are less expensive than human clinical trials, canine oncology clinical trials are significantly more costly than pre-clinical trials in mice [20]. There are also logistical difficulties in collecting 'fresh' versus formalin-fixed canine tumour specimens. Both these factors continue to limit the numbers of canine tumours subjected to genome-wide molecular profiling. However, it is the identification and characterisation, in canine tumours, of potential human tumour metastasis drug targets that is the basis for testing anti-metastasis agents on canine cancer patients. Thus, more appreciation of the clinical benefits of comparative oncology research and the use of companion animals in oncology clinical trials would lead to more financial support for such trials by government, philanthropic organisations and pharmaceutical companies. The recognition that there are advantages to using a One Health approach to developing new anti-metastatic therapies may lead to opportunities to fund such studies by cutting across the entrenched barriers of funding human-only or dog-only research studies and clinical trials and allow integrated projects to be undertaken.
A major hurdle for the development of anti-metastatics in either species is the Response Evaluation Criteria In Solid Tumours (RECIST), which is a set of published guidelines, developed in 2000 and updated in 2009 (RECIST 1.1), that define whether tumours in cancer patients shrink, stay the same, or enlarge during treatment [28]. Applying these response guidelines, the success of a drug in clinical trials relies on the demonstration of tumour shrinkage by X-rays, CT scans or MRI scans, with a confirmatory improvement in clinical outcome, and as such does not consider the capability to inhibit the development of metastasis as a measure of success. Only success in shrinking existing primary or metastatic lesions allows a drug to proceed to clinical trials for use in the adjuvant setting, with the aim of preventing/delaying the development of new metastases. Cancer Research UK, Cancer Research Technology and Cancer Therapeutics CRC Australia formed a Metastasis Working Group with representatives from academia, industry, government and regulatory bodies to develop recommendations on how to tackle the challenges associated with treating metastatic disease and reported that "…successful development of effective anti-metastatic therapies will require the regulatory agencies to work together with researchers, drug developers and statisticians to redefine the clinical development paradigm in order to encourage development of this complex but high-potential category of oncology drugs" [1]. Clinical trials in companion canines could help in this respect by incorporating the prevention of the development of metastases as a measurable clinical outcome that is given consideration when calculating the impact of a therapy. Considering that more than half of all dogs with appendicular OSA will have metastases by one year after diagnosis, the timeline for obtaining meaningful outcomes in an anti-metastatic canine trial is much shorter than in humans. It is hoped that pharmaceutical companies and regulatory bodies would then see the benefits that this offers and consider either (i) repurposing drugs classified as 'failed' by RECIST measurements in anti-metastatics clinical trials and/or (ii) change the RECIST measurements to include consideration of the inhibition of development of metastases as an outcome of importance.
Another key hurdle that the development of anti-metastatics must overcome is their high failure rate in human clinical trials, making anti-metastasis drug development appear to be a poisoned chalice. One possible reason for this high rate of failure is the fact that pharmaceutical companies do not routinely test drug candidates in metastatic models before moving them into clinical trials. Therefore, utilising a pathway to human clinical trials that moves from 'preclinical research' on cell lines and mouse models to 'veterinary clinical trials', utilising companion animals as predictors of human efficacy studies, has been proposed [29]. Such an approach opens the door for companion canines with a spontaneously occurring primary tumour and metastases, or no detectable metastases, to be involved in clinical trials of anti-metastatic drugs. The results of such trials may allow better prioritization of which drugs to take to human trials, and could thus potentially result in fewer failures. With more 'home runs' on the board, the prospect of increasing successful anti-metastasis drug development will be more appealing.

Conclusion
The similarities between the aetiology, pathogenesis and clinicopathological characteristics of human and canine cancers afford a rationale for both cross-species research to expedite the identification of targets for anti-metastatic agents, and the pre-clinical testing of anti-metastatics on canine cancer patients. The generosity and enthusiasm of dog owners to enable comparative research, either through sample donation or engagement in veterinary clinical trials, is a great strength of this strategy. The translational success of the companion canine metastatic OSA model has demonstrated the efficacy of veterinary clinical trials in facilitating the development of anti-metastasis drugs for human patients. Given the paucity of success in human clinical trials of anti-metastatics tested in conventional preclinical tumour models, companion canines with metastatic cancers represent a unique preclinical 'model' that is significantly under-utilised.
Open Access This article is distributed under the terms of the Creative Commons Attribution 4.0 International License (http://creat iveco mmons .org/licen ses/by/4.0/), which permits unrestricted use, distribution, and reproduction in any medium, provided you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons license, and indicate if changes were made.


extfieldsofstudy: Biology Medicine

text:
Taxonomic integrative and phylogenetic identification of the first recorded Triatoma rubrofasciata in Zhangzhou, Fujian Province and Maoming, Guangdong Province, China

Background Most species of Triatominae live exclusively in Latin America. However, one species, Triatoma rubrofasciata, has been recorded in the Americas as well as in various port areas in Africa and Asia. An increasing number of T. rubrofasciata have been reported in southern China in recent years. However, the origin of this invasive insect vector in China remains unknown, therefore, accurate identification and phylogenetic analysis of the bugs are urgently needed. Methods A total of seven triatomine insect specimens were found and collected from Maoming City, Guangdong Province, China (GDMM) and Zhangzhou City, Fujian Province, China (FJZZ), respectively. The obtained insect vector specimens were observed under a dissecting microscope for morphological classification and then the genomic DNA was extracted, and the 16S ribosomal RNA (rRNA), 28S rRNA as well as cytochrome oxidase subunit I (COI) genes of the species were amplified and sequenced. Subsequently, molecular phylogenetic analyses based on multiple alignments of the above genes were conducted in order to identify the species and determine the phylogenetic origin approximation accurately. Results The triatomine insects collected from GDMM and FJZZ were identified as Triatoma rubrofasciata using morphological and genetic analyses. All of the Chinese T. rubrofasciata captured in FJZZ, GDMM and other localities in southern China, together with a Vietnamese and Brazilian strain, formed a new, cohesive clade. T. rubrofasciata in GDMM and FJZZ are likely derived from strains found in Vietnam or Brazil. Conclusions To the best of our knowledge, this is the first record of the invasive insect T. rubrofasciata, which is likely derived from strains native to Vietnam or Brazil, in both Maoming City, Guangdong Province and Zhangzhou City, Fujian Province of China. A comparison of the DNA sequences of the 16 s rRNA, 28 s rRNA and COI genes confirmed the specific identification of T. rubrofasciata, and its potential origin in China is based on the phylogenetic analyses undertaken in this study. More targeted interventions and improved entomological surveillance are urgently needed to control the spread of this haematophagous insect in China. Electronic supplementary material The online version of this article (10.1186/s40249-019-0579-8) contains supplementary material, which is available to authorized users.

Multilingual abstracts
Please see Additional file 1 for translations of the abstract into the five official working languages of the United Nations.

Background
Triatomines, also known as kissing bugs due to their tendency to bite human faces, are insects from the order Hemiptera and family Reduviidae; they are the main vectors responsible for the transmission of American trypanosomiasis (Chagas disease), a chronic, systemic, tropical parasitic disease caused by Trypanosoma cruzi that infects an estimated 8 million people worldwide, mostly in Latin America [1][2][3]. There are a total of 151 triatomine species (149 extant and two extinct) that are currently composed of five tribes: Aberproseniini, Bolboderini, Cavernicolini, Rhodniini and Triatomini. Most Triatomine occur in the New World; only one genus, Triatoma, which is included in the Triatomini, is found in both in the Old and New Worlds [3,4].
Triatoma rubrofasciata is characterized by a triangular scutellum and an orange-red margin along the outer edge of the abdomen and the sides of the pronotum [5]; it is a cosmopolitan species from the genus Triatoma and is widely distributed across the globe, although it is most commonly found in Asia, Oceania, Africa and Central America [6][7][8]. This tropicopolitan insect is also widely spread across the southern coastal areas of China including Guangdong, Guangxi, Hainan and Taiwan [6,9]. Although Chagas disease has not been reported in China, dermatitis and anaphylactic shock caused by triatomine insects have been recorded since the 1980s [10,11]. Trends in global vector migrations from endemic to non-endemic countries have greatly increased the risk of Chagas disease expansion to global levels [12]; migration of these insect vectors could help spread this disease beyond its previous limits [13], which suggests there is a potential risk of Chagas disease being introduced and transmitted in China, posing a significant threat to public health [14]. Therefore, accurate identification and phylogenetic analysis of T. rubrofasciata will trace back its potential origin and carry out scientific and effective prevention, which may helpful for the control of its spread.
The distinctive features of the triatomine are based primarily on morphological characteristics such as the genitalia [15,16] and appendages [17,18], although molecular markers have also been used. However, morphological identification is complicated, time-consuming and sometimes inaccurate, especially for specimens in nymphal stages, due to the wide variety of triatomines (over 140 species). Recently, with the rapid development of molecular biology techniques and bioinformatics, molecular identification based on specific DNA sequence analyses has become a common technological tool used to identify new species [19,20], to provide rapid diagnoses of invasive pathogens [21] and cases of rare diseases [22], and to help in understanding the phylogenetic and phylogeographic relationships among species [23]. DNA markers including mitochondrial 16S and nuclear 28S ribosomal RNA (rRNA), mitochondrial cytochrome oxidase subunit I (COI), as well as DNA markers of nuclear origin, such as the internal transcribed spacer 2 (ITS-2) have been used to analyse genetic variation among triatomines, and explain their population structure and evolutionary history [6,10,[23][24][25][26][27][28]. Therefore, this paper aims to use the above three universal markers for the identification of triatomine species found in Maoming, Guangdong Province and Zhangzhou, Fujian Province, China and to trace their potential origin through phylogenetic analysis.  (Fig. 1). The locations chosen for collection were based on oral reports from local residents who found the insects in their houses. Triatomines were captured and transported live to a laboratory and identified using morphological features under a stereoscopic microscope (SZ650, Cnoptec, Chongqing, China) as described previously [29].

DNA extraction, amplification, cloning and sequencing
Total genomic DNA was isolated from three legs of each triatomine specimen using the HiPure Tissue DNA Micro Kit (Cat No. D3125-02, Magen, Guangzhou, China) according to the manufacturer's recommended protocol. The genes of mitochondrial 16S rRNA, COI and nuclear 28S rRNA were amplified using primers synthesized by Sangon Biotech (Shanghai, China), which are listed in Table 1.
Polymerase chain reactions (PCR) were performed in 50 μL of the following mixture: approximately 30 ng of DNA template, 1.1 × Golden Star T6 Super PCR Mix (Tsingke, China) and 0.4 μmol/L of each primer using a Bio-Rad PCR C1000 Touch instrument (Bio-Rad, USA). Negative controls were made with ultrapure water to control for contamination during PCR. For the 16S rRNA, 28S rRNA and COI genes, the fragments were amplified using the following thermal cycling conditions: initial denaturation at 98°C for 2 min; 35 cycles of denaturation at 98°C for 10 s, annealing at 55°C for 30 s, extension at 72°C for 15 s and final extension at 72°C for 3 min. Afterwards, the PCR products were separated using 1% agarose gel electrophoresis, and the expected bands were excised under ultraviolet light, followed by purification with an Agarose Gel DNA Extraction Kit (Takara, Japan). The purified DNA fragments were inserted into a pClone007 Blunt Simple Vector (Tsingke, China), and the plasmids were used to transform E. coli DH5ɑ. Sequencing was performed using the vector primers (M13 forward and M13 reverse) in an Automated

Sequence alignment and phylogenetic analysis
The sequences of the 16S rRNA, 28S rRNA and COI genes were submitted to the National Centre for Biotechnology Information (NCBI) GenBank database (https://www.ncbi.nlm.nih.gov/genbank/) under the accession numbers MK601647, MK601646, MK602653, MK602652, MK614012 and MK614011, respectively. They were submitted to GenBank for homology searches subsequent to triatomine species classification using BLAST (http://blast.ncbi.nlm.nih.gov/Blast. cgi) [30]. For the three genes, multiple alignments were performed among the best hit sequences using ClustalW2 [31] with default parameters; phylogenetic trees were constructed with MEGA 7.0 [32] using the Maximum Likelihood (ML) method [10] based on the Kimura 2-Parameter (K2P) nucleotide substitution model [33] and assessed with the bootstrap-resampling technique over 1000 replications. To demonstrate the polymorphism of T. rubrofasciata from multiple localities and their intraspecific genetic divergence, the pairwise genetic distances of 16S rRNA, 28S rRNA and COI were calculated between sequences obtained in this study and reference sequences downloaded from GenBank using the K2P nucleotide substitution model in MEGA 7.0 (www. megasoftware.net) [32].

Morphological characteristics of Triatoma rubrofasciata
Four adult triatomine insects (two males and two females) were caught on a balcony where lumber was stored on the second floor of two farmers' dwellings at Huping Village, Fengshan Town, Hua'an County (FJZZ). In addition, three adult triatomine insects (two males and one female) were captured under wooden beds in a living-room on the second floor of a local resident's house at Tiantou Village, Shuidong Town, Dianbai District (GDMM). The triatomines collected in FJZZ and GDMM showed the morphological characteristics of T. rubrofasciata (Figs. 2 and 3). The specimens were black with orangered margins along the lateral edge of the pronotum, wings and abdomen. The forepart of the pronotum was narrow, whereas the posterior was wide. The female and the adult males were approximately 23 mm and 21 mm in length respectively, with slightly elliptical shapes, and the long, cone-shaped heads were approximately three millimeter in length (Fig. 2). The thread-like antennae, located between the black compound eyes and clypeus, consisted of four segments, with the 1st segment surpassing the apex of head (Fig. 3a, b). Numerous small particles were densely distributed around the head and pronotum (Fig. 3a). The scutellum was triangular in shape with wrinkles in the middle and a sharp, needlelike apex (Fig. 3a). The genital segment of the adult female was triangular form, while that of the adult male had a broad circular shape (Fig. 3c, d).

Molecular identification and phylogenetic relationships of Triatoma rubrofasciata
In this study, apart from morphological identification, we also used a molecular identification method based on sequences of 16S and 28S rRNA as well as COI to classify the triatomine species from GDMM and FJZZ. The results of the similarity analyses on the basis of the above genes demonstrated that the triatomine species in these two localities were accurately identified as T. rubrofasciata with 96-100% sequence matching identities.
For the 16S rRNA gene, T. rubrofasciata from FJZZ displayed the highest similarity with specimens from Foshan City, Guangdong Province, China (GDFS) with 99% identity. T. rubrofasciata from GDMM demonstrated the highest similarity with T. rubrofasciata from Zhanjiang City, Guangdong Province, China (GDZJ) with 99% identity. A polymorphism analysis of the 16S rRNA gene showed that two transition sites differed between the sequences of the insects from FJZZ and GDMM, and from FJZZ and GDZJ. Similarly, two distinctive transition sites were observed between the 16S rRNA sequences of the insects from GDMM and GDFS. Furthermore, one transversion site and one transition site differed between the 16S rRNA sequences of T. rubrofasciata from FJZZ and TW, while one transversion site and three transition sites differed between the sequences of the insects from GDMM and TW. However, only one distinguishing transition site was detected between the sequences of the insects from FJZZ and Vietnam and from GDMM and Vietnam. The genetic distance of the 16S rRNA gene ranged from 0.000 to 0.003 in T. rubrofasciata from the different locations, and the average intraspecific distance was 0.001 (Table 2, Additional file 2: Table S1).
For the 28S rRNA gene, both T. rubrofasciata from FJZZ and GDMM were similar, with T. rubrofasciata from GDZJ demonstrating 100% identity. The 28S rRNA sequences of T. rubrofasciata from FJZZ and GDMM were identical and only one distinctive transversion site was found between the sequences of T. rubrofasciata from FJZZ and GDFS, whereas no sites differed between the sequences of T. rubrofasciata from FJZZ and Ningming County, Guangxi Province, China (GXNM).
In addition, two different transversion sites and one distinguishing transition site were found between the sequences of T. rubrofasciata from FJZZ and Brazil, while two transversion sites differed between the sequences of the triatomine insects from FJZZ and Vietnam. The genetic distance of the 28S rRNA gene in T. rubrofasciata from the different locations ranged from 0.000 to 0.002, and the average intraspecific distance was 0.001 (Table 3, Additional file 3: Table S2).
The COI gene of T. rubrofasciata from both FJZZ and GDMM presented the highest similarity with that from Brazil with 96% identity, with one variable site between them. Nevertheless, eighteen variable sites were detected between the COI sequences of the insects from FJZZ and Brazil, while nineteen variable sites were observed between the sequences of the insects from GDMM and Brazil. Thus, the genetic distance of the COI gene varied from 0.003 to 0.048, and the mean intraspecific distance was 0.032 (Table 4, Additional file 4: Table S3). Furthermore, phylogenetic analysis of triatomine species based on the sequences obtained in this study and reference sequences with high similarity obtained from GenBank for the 16S rRNA, 28S rRNA and COI genes showed a distinct clade for all T. rubrofasciata individuals and was strongly supported by bootstrap values as high as 97 to 100%, which also confirmed that the insects captured in FJZZ and GDMM were T. rubrofasciata (Figs. 4, 5 and 6). Moreover, these results indicate that the triatomine species in the ML phylogenetic trees might be identified accurately using the barcode region used in this study.
According to the phylogenetic trees inferred from the 16S rRNA and 28S rRNA genes (Figs. 4 and 5), the Triatomini tribe is scattered throughout 4 main clades: (1) Clade A: South American Triatoma, including the brasiliensis subcomplex (T. brasiliensis from Brazil), the rubrovaria subcomplex (T. rubrovaria and T. circummaculata from Argentina), the infestans subcomplex (T. infestans from Argentina) and the matogrossensis subcomplex (T. matogrossensis and T. costalimai from Brazil); (2) Clade B: genus Panstrongylus (P. geniculatus and P. megistus from Argentina, and P. lignarius from Brazil); (3) Clade C: North and Central American Triatoma, including the protracta complex, (T. protracta from the United  Fig. 5. This strongly suggests a close relationship between Asian triatomines and the rubrofasciata complex from Brazil. Similarly, phylogenetic analysis of Triatominae COI sequences demonstrated that Chinese triatomines from FJZZ and GDMM and the Brazilian T. rubrofasciata were in the same clade (Fig. 6). The Triatomini tribe was scattered throughout three clades of the South American Triatoma (genus Panstrongylus) as well as the South American and Asian Triatoma, due to the limited sequences of the COI gene available in GenBank (Fig. 6). These results strongly imply that Vietnam or Brazil is the most likely origin of the T. rubrofasciata in China.

Discussion
T. rubrofasciata is a widespread haematophagous insect [34,35], and it is the vector for T. cruzi, which causes Chagas disease and is a potential vector of bacteria and viruses [36]. It is likely distributed via shipping lines from the tropics and sub-tropics [34,37,38], and from inland regions in north-east India, Vietnam [35] and north-east Brazil [39]. Though T. rubrofasciata is the first species of the triatomine to have been described, its origins remain a mystery [40]. Several researchers have assumed an Old World (oriental) origin for T. rubrofasciata [29,41]. However, more recent studies have suggested that all Triatominae are of a New World origin and that they were transported in association with rats (Rattus norvegicus) on ships sailing from the Americas [42][43][44]. Clear identification of invasive species is crucial for ecosystem maintenance, control, and risk assessment [19]. In this paper, the triatomine were identified as T. rubrofasciata by morphological analysis, which is supported by both the genetic and phylogenetic species concepts, indicating that 16S rRNA, 28S rRNA and COI can be used as specific gene markers for the reliable identification of triatomine species, as suggested by previous studies [6,10,44]. Morphological taxonomy and molecular identification are the two primary methods used for the identification of species [45][46][47]. Morphological similarities at both the intra-and inter-species level limit the usefulness of morphological taxonomic keys in some species [48]. Moreover, molecular identification might be effectively applied to classify the immature stages of triatomines, as traditional taxonomic methods are difficult to use when distinguishing nymphs, especially considering the species complexes that exhibit large morphological and chromatic variation [28]. However, molecular identification will not take the place of morphological systematics because it relies on species definitions produced from morphological analyses [49]; hence, both approaches should be applied to provide effective, efficient and accurate species identifications [50,51], and the integration of morphological and molecular methods is indicated for species recognition of T. rubrofasciata in future entomological surveys in China [6].
To determine genetic variation in T. rubrofasciata from the different sites, genetic distances were calculated. The polymorphism analyses results showed that transitions  were more frequent than transversions in the 16S rRNA, 28S rRNA and COI genes, and our data was consistent with observations from previous studies suggesting that COI is a faster evolving sequence, exhibiting higher variability than 16S or 28S rRNA [52,53]. Hence, we suggest that 16S rRNA, 28S rRNA as well as COI genes are suitable for species identification while for tracking back the origin of the species, COI gene is applicable due to the higher genetic diversity. Despite a relevant body of research on taxonomy, distribution, pathology and prevention [13], little is known about the phylogeny of T. rubrofasciata in Asia, as DNA sequence information regarding this vector is scarce. While an increasing number of T. rubrofasciata have occurred in southern China [6,[9][10][11][12], their origins remain unclear. Maximum Likelihood phylogenetic trees inferred by 16S rRNA and 28S rRNA gene sequences suggest that the Triatoma could be divided into four clades, not three as suggested by a previous study that used a limited sample size [9], and that the T. rubrofasciata from China, Vietnam and Brazil form a new, cohesive clade, indicating a close relationship between T. rubrofasciata from Asia and South America, with no incipient speciation between them.
Though the interrelationships among genera and species of the Triatomini have been shown as unstable in previous studies [43,44], most researchers consider the Asian species of T. rubrofasciata to be derived from the Americas and form a monophyletic clade within the Triatominae [4,44], which is consistent with the findings from this study. The results of the comparative analyses of the 16S rRNA, 28rRNA and COI sequences of T. rubrofasciata exhibited high genetic similarity, suggesting a very recent dispersal of this invasive insect in Asia (including China and Vietnam) and a derivation of the insect from South America, primarily from Brazil.
Though this study does allow a better understanding of the phylogenetic relationships among triatomine species and determine the phylogenetic origin approximation of T. rubrofasciata in China, extensive   T. rubrofasciata identified in this study were shown in the square frame population genetic study is needed to confirm these results due to the small sample size and limited sampling area in this paper. Especially for the COI gene, we can only obtain the sequences of T. rubrofasciata from FJZZ and GDMM in China as well as Brazil, but not Vietnam, from NCBI GenBank database. Moreover, the morphological difference between the bugs from China, Vietnam or Brazil is missed in the present study because the authors have not the samples of T. rubrofasciata captured in Vietnam and Brazil.

Conclusions
This is the first report of the triatomine species found in GDMM and FJZZ in China, which likely originated in Vietnam or Brazil. The current molecular taxonomic study and phylogenetic origin analysis revealed genetic data of T. rubrofasciata from different localities in China, which could contribute to the establishment of adequate vector control strategies and adding phylogenetic information regarding this widespread, invasive insect of medical significance.

Additional files
Additional file 1: Multilingual abstracts in the five official working languages of the United Nations. (PDF 597 kb) Additional file 2: Table S1. Genetic distance table for 16S rRNA sequences used in Fig. 4. Standard error estimates were shown above the diagonal. Analyses were conducted using the Kimura 2-Parameter model. (XLSX 26 kb) Additional file 3: Table S2. Genetic distance table for 28S rRNA sequences used in Fig. 5. Standard error estimates were shown above the diagonal. Analyses were conducted using the Kimura 2-Parameter model. (XLSX 20 kb) Additional file 4: Table S3. Genetic distance table for COI sequences used in Fig. 6. Standard error estimates were shown above the diagonal. Analyses were conducted using the Kimura


s2fieldsofstudy: Biology Medicine

text:
Identification of Targetable Pathways in Oral Cancer Patients via Random Forest and Chemical Informatics

Treatment of head and neck cancer has been slow to change with epidermal growth factor receptor (EGFR) inhibitors, PD1 inhibitors, and taxane-/plant-alkaloid-derived chemotherapies being the only therapies approved by the U.S. Food and Drug Administration (FDA) in the last 10 years for the treatment of head and neck cancers. Head and neck cancer is a relatively rare cancer compared to breast or lung cancers. However, it is possible that existing therapies for more common solid tumors or for the treatment of other diseases could also prove effective against oral cancers. Many therapies have molecular targets that could be appropriate in oral cancer as well as the cancer in which the drug gained initial FDA approval. Also, there may be targets in oral cancer for which existing FDA-approved drugs could be applied. This study describes informatics methods that use machine learning to identify influential gene targets in patients receiving platinum-based chemotherapy, non-platinum-based chemotherapy, and genes influential in both groups of patients. This analysis yielded 6 small molecules that had a high Tanimoto similarity (>50%) to ligands binding genes shown to be highly influential in determining treatment response in oral cancer patients. In addition to influencing treatment response, these genes were also found to act as gene hubs connected to more than 100 other genes in pathways enriched with genes determined to be influential in treatment response by a random forest classifier with 20 000 trees trying 320 variables at each tree node. This analysis validates the use of multiple informatics methods to identify small molecules that have a greater likelihood of efficacy in a given cancer of interest.

Introduction
Head and neck cancers can be defined as cancers in the upper airway and/or digestive tract found in oral cavity and laryngeal, pharyngeal, oropharyngeal, and hypo-pharyngeal tissues. Such cancers make up 3% of cancers diagnosed each year. 1 Head and neck cancer incidence has declined from 25 cases per 100 000 in the 1990s to 15 cases per 100 000 in the present day. 2 While the decrease in head and neck cancer incidence may be due to a drop in tobacco use, 3,4 the mortality associated with head and neck cancers has not changed significantly in the last 20 years. Treatment of head and neck cancer has also been slow to change with epidermal growth factor receptor (EGFR) inhibitors, PD1 inhibitors, and taxane-/plant-alkaloid-derived chemotherapies being the only therapies approved by the U.S. Food and Drug Administration (FDA) in the last 10 years for the treatment of head and neck cancers. Head and neck cancer is a relatively rare cancer compared to breast or lung cancers. However, it is possible that existing therapies for more common solid tumors could also prove effective in oral cancers. Many therapies have molecular targets that could be appropriate in oral cancer as well as the cancer in which the drug gained initial FDA approval. Also, there may be targets in oral cancer for which existing FDA-approved drugs could be applied. This study describes informatics methods that use machine learning to identify influential gene targets in patients receiving platinum-based chemotherapy, non-platinum-based chemotherapy, and genes influential in both groups of patients.
Drugs approved by the FDA for oral cancer are methotrexate, cetuximab, pembrolizumab, nivolumab, and docetaxel. These therapies are combined to be used in conjunction with platinum-based chemotherapies such as cisplatin or carboplatin unless those therapies are contraindicated due to comorbidities such as renal disease. The small number of new oral cancer drugs could be attributed in part to the low overall burden of oral cancer in comparison with other cancers. The current timeline for FDA approval of a novel small molecule or biologic is 10 years or more. Repurposing existing FDA-approved drugs is a popular method used to shorten this process to 3 to 4 years. Cetuximab is an EGFR inhibitor that has been shown to decrease the rate of progression of oral cancer used in conjunction with Cisplatin. Pembrolizumab and nivolumab are both PD1 inhibitors that used T-cells to attack cancer while it progresses. These therapies use the body's immune system as another treatment modality to reduce the burden of oral cancer. Current literature provides support for the role of ligand channel gating, hedgehog signaling, 5-7 NOTCH, B-WICH, 8 inflammasome, 9 WNT, 10 and Calci neurin pathways in cancer. [11][12][13] The role of these pathways and targeting specific genes within them has been pursued in other cancers, but have with few exceptions, not yet been examined in oral cancer. Possible gains from targeting these pathways would be initiating immune response, targeting cancer metabolism, targeting signaling for metastasis, and targeting inflammation pathways that may drive progression. If there is a synergistic effect to attacking 2 Cancer Informatics multiple hallmarks of cancer simultaneously, then the net gain to the patient would be a net gain in overall post diagnosis survival time for the oral cancer patient.
Identifying a means by which drugs may be prioritized for further screening and validation for a specific cancer type would be desirable. Databases linking genes, the proteins they express, ligands corresponding to those proteins, and structural data that can be analyzed all exist in varying forms or completeness across different publicly available databases. This study describes how integration of analyses of these databases can be used to select gene targets in a specific cancer and how therapies can be prioritized for screening based on existing structural information for the ligands associated with genes and the proteins they express.
There are several hurdles to the analysis of high-dimensional genomic data using traditional regression analyses. Random forest analysis is a machine learning approach that is less hindered by datasets with large predictor to observation ratios. In this study, we apply random forests analysis to gene expression data to identify those genes and pathways that are most predictive of post diagnosis survival across treatment strata. This is the first application of this approach to head and neck cancer patient data in The Cancer Genome Atlas (TCGA). National Comprehensive Cancer Network (NCCN) 14 guidelines recommend that node-positive patients with tumors of clinical stage 3 and greater receive chemotherapy. Platinum-based chemotherapy with radiation and surgery is the current standard of care recommended to these patients. Patients that do not receive chemotherapy recommendations by NCCN are node negative with clinical stage 2 and lower. Following standards set by NCCN guidance on treatment, this study chose to identify patients receiving and not receiving platinum-based chemotherapy as separate groups. Analysis of influential genes in each group will improve knowledge of possible mechanisms driving treatment response for early stage and more advanced tumors.
Random forest is a machine learning approach to identifying the most important predictors in high dimensional datasets.15,16 This approach is uniquely suited for classification of observations in datasets where P (predictors) are greater than N (number of observations). Random forest randomly selects predictors from a large group of predictors and then applies those predictors to a decision tree predicting overall survival. Random forest does not pay a statistical penalty when the number of observations is small. Instead the strength and limitation of this method is its reliance on computational intensity. That is, as the number of decision trees in a random forest increase, so does classification accuracy. Accuracy is also dependent on the number of predictors tried at decision tree nodes. As node size and forest size increase, so does forest classification accuracy. However, there is a rate of diminishing returns in the accuracy gained from each tree added to a forest. Therefore, computational time and cost must be factored into all random forest analysis plans to measure project feasibility. Random forest has been successfully applied to predicting cancer diagnosis and treatment response for a variety of cancers. [17][18][19][20][21] For this study, we have selected to apply random forest analysis to the gene expression values of oral cavity cancer patients to identify the upregulated pathways most predictive of improved treatment response across gender and environmental exposure subgroups like alcohol and tobacco. RNAseq data are inherently high dimensional, applying typical regression models to such data can be costly as large sample sizes are required to identify even moderate effect. Identifying gene interactions can be even more costly in terms of the required statistical power. Stratified pathway analysis via random forest methods has been shown to be successful in identifying single influential genes (within the context of larger pathways) that are predictive of overall survival with limited sample size. 22 This approach has not yet been applied to identification of influential genes and gene interactions within oral cancer patients stratified specifically by treatment. In this way, the importance of pathways and genes of interest can be compared across strata to assess which subgroups may be most sensitized to changes in gene expression within a given pathway.

Methods
This study focuses on the identification of the role of gene expression in oral cavity cancer patients and applying machine learning approaches like random forest to determine genes that are important in influencing treatment response. Reference ligands known to bind to proteins expressed by genes deemed influential by random forest can be sent through a virtual screening pipeline to identify small molecules with greater likelihood of acting as protein agonists/antagonists. Ligands that have a strong shape similarity to known binding ligands have greater potential for success in high-throughput screening endeavors. As shape similarity alone is insufficient in identifying new drug leads, all leads will be validated with existing literature, and those leads without previous biological validation will be presented as such.
By using a stratified random forest analysis, we will be able to rank genes within the strata of chemotherapy treatment status. This approach will allow for the identification of those top ranked genes that are unique to each stratum. This will be done by identifying common and unique genes between sets of genes influencing the treatment response in patients receiving platinum-based chemotherapy and those that do not. The result will be the identification of oral cavity cancer pathways influencing treatment response which will inform researchers on mechanisms driving treatment response in specific groups such as late-stage, node-positive patients who are more likely to receive chemotherapy treatment. This analysis will illustrate and support existing studies showing the strength of machine learning methods as an alternative method in identifying gene expression values influencing treatment response. This study is focused not only on the predictive power of an aggregated panel of gene expression values, but rather to integrate random forest with chemical informatics and thus describe methods to shorten the pursuit of novel therapies treating cancers with relatively lower incidence.

Retrieval of Public Data
This study used clinical and genetic data obtained from TCGA. Genetic data included raw counts per million (CPM) of RNA sequence expression values for 523 patients posted to TCGA. Of these 523 patients, 313 were diagnosed with Oral Squamous Cell Carcinoma (OSCC). Oral squamous cell carcinoma patients included tongue, buccal mucosa, alveolar ridge, general oral cavity, and soft palate tissues. Of these 313 patients, 267 were included based on complete survival time. Of these 267, 109 received either Carboplatin, Oxaliplatin, or Cisplatin, while 158 patients received a treatment other than the platinum-based chemotherapy treatment regimen. All tissue samples were collected prior to start of treatment. Clinical data on tumor stage, necrosis, size, and nuclei were also retrieved from TCGA. We obtained demographic data on ethnicity, race, and gender from TCGA clinical files. In addition, this dataset had information on environmental exposures like tobacco history (ever/never smoke) and number of alcohol containing drinks in a day (greater than 2 drinks consumed per day, 2 drinks or less consumed per day). Overall survival time in months was extracted as a measure of treatment response.

Description of random forest approach
Stratified pathway analysis considers important covariates in data analysis. In this analysis of head and neck cancer gene expression data, this study used information regarding patient age, sex, alcohol, and smoking exposures. In the first round of analysis, the sample of 267 OSCC patients was organized into 109 and 158 subsamples based on whether patients received platinum or nonplatinum-based chemotherapy, respectively. For each group, this study built a random forest 15 to predict survival time based on the gene expression levels within each subgroup.
To better communicate the function of random forests, understanding of a decision tree construction is needed. A decision tree is constructed by the following steps: Step A: A bootstrap sample was taken from the original sample.
Step B: A decision tree is grown for each bootstrap sample.
Step C: At each tree node, a predetermined number of predictors were applied randomly to create branches within the tree.
Step D: A branch is formed using the predictor from step C.
Step E: Steps C and D are repeated until the end of every tree branch contains samples above or below the same survival threshold or contain only one sample.

Random forest result measurements
Random forests build many decision trees to comprise a forest. Each tree is put together by using a random bootstrap sample of the original data and applying a random number of predictors at each node of a decision tree. The SRCRandomForest R package 23 employed in this study sets aside half of the data to be used for validation purposes to measure the accuracy of the random forest model. The SRCRandomForest package was chosen for this study because it applies a multivariate Cox regression model to produce each decision tree within each Forest. In this way, censored and non-censored data can be considered when carrying out this analysis. In this way, the model can measure the degree of influence of each gene on patient survival. The P values yielded through this analysis are defined as the "proportion of cross-validation errors smaller than the cross-validation errors obtained from 500 iterations of random forest runs of randomly permuted labels of patients." This list of genes can be used to identify pathways that are enriched with the influential genes identified through random forest at odds that would be greater than can be attributed to chance alone with a P value of .05. This analysis will present pathways common to, and unique to, each chemotherapy treatment strata. Such analyses may identify plausible biological mechanisms that enhance understanding of observed differences in survival. To reiterate, the focus of this study is not to pursue a diagnostic tool but to identify those gene expression values exerting a strong influence on treatment response. This study adopted random forest as the machine learning method of choice due to its superior interpretability and scalability. Random forest and many applications incorporating random forest have been developed more than 30 years. Random forest was of interest to this study as an application developed specifically for the use of survival data had been developed as an open-source R package. Although random forest has been shown to be less efficient than extreme gradient boosting applications, the long history of this approach has allowed it to evolve to produce an application well suited toward achieving the analytic goals of this study.

Random forest tuning parameters
Our random forest model used 20 000 tree forests for a forest size, with 320 variables tried at each node in each decision tree; 20 000 trees were the point at which we could identify no significant increase in our ability to predict patient survival; 320 tries at each node in each decision tree were double the recommended number of tries given by the author of our R software package RandomForestSRC. The author H.
Ishwaran et al 23 refers to a generally accepted practice of using the "square root of the total number of predictors as a starting point for the number of variables tried at each node." It was for this reason that we applied 20 000 trees and 320 variables tried at each node in each decision tree for every stratum in our analysis. This approach was applied to our entire final sample of 267 OSCC patients. We then divided this sample by whether a patient received platinum-based chemotherapy or not. All random forest analyses were carried out on 64 computational nodes in the University of California High Performance Computational Cluster. Random forest analysis was made more efficient by running analysis in conjunction with the R package parallel. The output of each group's analysis produces a list ranking each gene. This analysis identified common and unique pathways between the entire dataset and each chemotherapy treatment group. We then identified the unique and common genes between chemotherapy groups. This analysis will allow us to observe the difference in gene importance and corresponding pathways in relation to overall survival.

Description of virtual screening approaches
Using chemical informatics techniques, ZINC drug database 24 of 1379 FDA-approved drugs (FDA) and ZINC Traditional Chinese Medicine (TCM) database of 39 894 small molecules can be used to apply three-dimensional chemical informatics approaches to the identification of small molecules that are the best candidates for inhibition of proteins expressed by those genes influencing treatment response. Reference ligands for each protein are obtained from the Royal Chemistry Society Protein Database 25 and then virtually screened against FDA and TCM small molecule libraries. Molecular shape overlay is an approach for measurement of the similarity of one molecule in comparison with another. A Tanimoto coefficient is used to measure the degree of similarity between 2 molecules.
The Tanimoto coefficient between 2 points, a and b, with k dimensions is calculated as follows The Tanimoto similarity only applies to binary variables, for binary variables the Tanimoto coefficient ranges from 0 to 1 (where 1 is the highest possible similarity). 26 A goal of this study is to perform searches of 2 small molecule databases FDA and TCM using a maximum common substructure measurement of Tanimoto similarity from the R Rcpi package 27 that has been shown to perform robustly across a variety of molecule types. Table 1 provides information on those genes that are significantly enriched within our gene set beyond what would be expected by chance alone. The significance of enrichment is calculated as the odds of randomly selecting the number of genes in the submitted set of genes by randomly selecting from 20 530 genes more than 100 times. The false discovery rate (FDR) reported by the Pathway Reactome application used for this study represents an adjustment for multiple comparisons across all pathways. Benjamini-Hochberg FDR is calculated as the P value ranking (smallest being 1 and all following having greater or equal rank dependent on size of P value) divided by the number of tests performed and multiplied by the significance criterion. In this study, .05 is the criterion used to measure significance. The FDR can be interpreted as the proportion of tests within a set of tests that falsely rejects the null hypothesis. If an FDR is .5, then 50% of those pathways identified will falsely reject the null hypothesis. It is important to note that the FDR calculation used by Pathway Reactome defaults to perform a large number of analyses/tests. In addition, there are many pathways examining similar genes and gene types. Unfortunately, this thorough examination strategy also inflates the number of analyses and causes the FDR to become overly conservative.

Pathway analysis
To identify those genes that are most likely to be connected to influential pathways, 2 filters were applied to gene selection. First, the gene had to be in the top 5% of influential genes identified via random forest. Second, if the gene was not within the top 5% of genes, it could still be included within the analysis if it was within the top 40% of genes and was known to be connected via past experimental studies supporting gene network connections found in the Cytoscape/Pathway Reactome plugin database. Gene topology can be accessed by uploading topology from a given gene in Cytoscape and merging it with that gene's corresponding random forest importance values. Gene topology refers to the connectedness of a gene (supported by experimental results referenced in cytoscape) to other genes within a Cytoscape network. Thus, the selection of genes included those ranked as the most important by random forest (top 5%) or of moderate importance and high topology (shown through past experiment or literature search to be connected to 50 genes or more). For each set of analyses, 1000 influential genes were selected and 100 high topology genes were selected. In this way, integration of Cytoscape network analysis with random forest results allowed for identification of pathways significantly enriched with the genes identified by the random forest model. For a complete overview of the analyses in this study please reference the online Supplemental Materials.

Top pathways
This section will describe pathways uniquely influential in patient response to platinum-based chemotherapy and influential in Schomberg 5 response to non-platinum-based therapy. Influential pathways shared by both platinum-based chemotherapy users and nonusers will also be presented. An influential pathway will be defined as a pathway that is significantly enriched with genes that were in the list of top 1000 (5%) of most influential genes yielded by random forest analysis for platinum-based chemotherapy users, nonplatinum-based chemotherapy users, and for those pathways enriched with genes shared in common in lists of top 1000 genes for platinum-based chemotherapy users and non-users. The top 5% of genes were conservatively selected to produce greater certainty of the link between highly ranked genes, the pathways in which they were enriched, and the link between treatment response and significant pathways identified through gene enrichment analysis. Top pathways could also be enriched with those genes in the top 40% of important genes identified by random forest if they also had a high amount of connectedness (a gene was connected to 25 gene nodes in a gene network) reported by the Cytoscape/Reactome application.
Pathways significantly enriched with genes identified as the most influential (top 5%) in predicting treatment response for platinum-based chemotherapy users were those related to calcium channel gating, hedgehog signaling, histone acetylation, elastic fiber production, tRNA acetylation, hexokinase deficiency, inhibition of adenylate cyclase, and CLEC7A inflammasome pathway. It has been reported that calcium channel gating has been associated with multiple cancers. [28][29][30] There have also been recent studies evaluating the benefit of targeting histone deacetylation pathways in oral cancer. [31][32][33] The hedgehog signaling pathway has also been shown to signal progression in other cancers "Hh signaling has been shown to regulate the selfrenewal of CSCs in breast, glioma and multiple myeloma, and more convincingly in the maintenance of chronic myelogenous leukemia (CML) stem cells." [34][35][36][37][38][39] All significant pathways for platinum-based chemotherapy users are in Table 1.
Significant pathways for patients not using platinum-based chemotherapy were those related to B-WICH complex, TP53 pathway, fibroblast growth factor receptor (FGFR) pathway, potassium channel gating, and RNA polymerase chain elongation pathways and their epigenetic regulation. TP53 and FGFR pathways represent the expression of canonical oncogenes which have been shown to be cancer drivers and associated with the production of all cancers. [40][41][42][43] The B-WICH complex has been found to be linked to maturation of invadopodium in breast cancer and has been suggested as both a biomarker and target for cancer invasiveness. 8,44 Potassium channel gating has also shown to be a potential target for head and neck cancers due to this pathways association with immune response and treatment response. [45][46][47][48][49] RNA polymerase chain elongation and its role in transcription is a logical contributor to cancer progression and differentiation; however, the lack of specificity makes this a difficult pathway to target specifically in cancer cells. All significant pathways for non-platinum-based chemotherapy users are in Table 2.
For those pathways enriched with genes shared by users and non-users of platinum-based chemotherapy, it was found that pathways related to g-protein beta folding, nuclear factor of activated T-cells (NFAT) activation, and repression of Wnt pathway genes were enriched with genes influential in treatment response in both groups of patients. Repression of Wnt pathway genes may be done through targeting the sonic hedgehog pathway as previously outlined or through more direct means which have been researched in multiple other cancers. [50][51][52] Nuclear factor of activated T-cells proteins have been found to be associated with cancer progression in blood and  [53][54][55] These pathways influencing treatment response in both users and non-users of platinum-based chemotherapy can be seen in Table 3.

Important genes and biological implications
A visualization of pathways overlapping between users and non-users of platinum-based chemotherapy highlights the importance of several genes in a way that random forest analysis alone could not. By visualizing the 4 common pathways, it becomes possible to identify not only highly influential genes but also those genes that have the highest degree of connectivity to influential genes. Using annotation built into Cytoscape, 56,57 we can also identify existing small molecules used in cancer therapy that are not yet commonly used in oral cancer, and we can also observe those genes previously found to be associated with oral cancer. Genes found to be influential in oral cancer for patients receiving platinum-based chemotherapy with existing literature supporting the targeting of these genes in cancer were INSR, BRAF, and PSMB7 which are targeted by Ceritinib (regorafenib and dabrafenib) and bortezomib, respectively.
These drugs are not currently FDA approved for treatment in oral cancer. Genes found to be influential in oral cancer for patients not receiving platinum-based chemotherapy with existing FDA-approved chemotherapy drugs targeting the products of said genes are FGFR3, EGFR, PRKAA2, CSNK2A1, INSR, MET, CAMK2A, PSMB5, and PSMB1. There are multiple chemotherapy drugs targeting these pathways with 14 different drugs targeting EGFR alone. It should be noted that EGFR is a gene pathway being targeted in current oral cancer treatment. Sarafenib Tosylate, Pazopanib Hydrochloride, and Vadetanib all target the FGFR pathway specifically. Sunitinib Malate is unique in that it has been found to act on 4 different genes that were found by random forest to be influential in treatment response: FGFR3, CDNK2A1, PRKAA2, and CAMK2A. Again, we see that Certinib acts on a gene that is influential in both patient treatment groups, and that gene is INSR. Bortezomib, Carfilzobib, and Ixazomib all act on PSMB5, which is an influential gene in both Platinum and non-platinum-based therapy. PSMB5 and PSMB1 are both found to be within the top 5% of influential genes in random forest analysis and are genes that are significantly enriched within the sonic hedgehog pathway. Chemotherapy drugs and their relationship to genes in common influential pathways between users and non-users of platinum-based chemotherapy are visualized in Figure 1 for platinum-based chemotherapy users and Figure 2 for non-platinum-based chemotherapy users.
In addition to analysis of the intersection of existing cancer drugs and genes deemed influential by random forest, this study also looked at the intersection between gene topology within a pathway and random forest influence. This analysis identified CTNNB1, PLCG2, SHC1, UBA52, UBB, UBC, and HDAC3 as genes that meet filters of belonging to one of the 4 common enriched pathways, being a gene that is one of the top 5% of influential genes listed by random forest analysis and being connected to more than 100 genes within the 4 interconnected  Hexagon shapes are genes. Dark red genes are of greater influence based on random forest analysis results (within top 5% of influential genes), white genes do not fall within the criteria of being in the top 5% of influential genes. All 4 common pathways enriched in patients receiving platinum and non-platinum therapy (described in Table 3) were merged into one pathway using Cytoscape. genes were clustered according to the Reactome Pathway Plugin available via Cytoscape.

8
Cancer Informatics pathways (Figures 1 and 2). CTNNB1 mutations have been found to be predictive of lung and other thoracic cancers, [58][59][60][61][62] and PLCG2 and calmodulin knockdown have been shown to induce paclitaxel sensitivity in cervical cancer tumors. This may prove of use to oral cancer patients, which may be assigned to paclitaxel or other taxol regimen. 63,64 SHC1 has been shown to be a regulator of EGFR function and thus a potential target for multiple cancer types where EGFR is a key driver. [65][66][67] Ubiquitin genes UBA52, UBB, and UBC have been shown to be associated with several cancers and research is currently being pursued in targeting ubiquitin ligases to improve treatment response. [68][69][70][71][72] Histone-deacetylase genes specifically HDAC3 is shown to be a hub to several genes that are influential in platinum-based chemotherapy response genes that have also been associated with metastatic invasion in breast and pancreatic cancer ( Figure  1). Inhibition of HDAC3 was shown to impact signaling to cancer stem cells. This gene has been shown to be a regulator of apoptosis control exerted by TP53. [73][74][75][76] These genes were not shown to currently have any antibody or small molecule therapies targeting their action. High-level topology (more than 100 connected genes) and high random forest importance ranking should provide impetus for further research into the targeting of gene action in oral cancer.

Chemical informatics analysis of drug targets and leads
For those genes meeting topology and random forest filters, the known ligands of proteins expressed by each gene were identified though the RCSB protein data bank. Structural files of ligands were downloaded as .sdf files and uploaded into the chemical informatics R package Rcpi. Once loaded, each ligand had to undergo virtual screening against all FDA-approved drugs to identify existing FDA-approved drugs that may prove efficacious as therapeutic agents. Only those molecules with a Tanimoto similarity score >50% were included in results. A TCM small molecule database was also used as biologicderived small molecules are known to provide better shape overlay when screened against other biologic small molecules.
In addition, the molecules in the TCM database have been Hexagon shapes are genes. Dark red genes are of greater influence based on random forest analysis results (within top 5% of influential genes), white genes do not fall within the criteria of being in the top 5% of influential genes. All 4 common pathways enriched in patients receiving platinum and non-platinum therapy (described in Table 3) were merged into one pathway using Cytoscape. genes were clustered according to the Reactome Pathway Plugin available via Cytoscape. loss of ring structure is indicative of differences in influence of genes between patients receiving platinum and non-platinum therapies.

9
shown to be generally safe in people by merit of its long historical use in human populations. For CTNNB1, several ligands were identified via RCSB PDB (2s)-3-{[{[(2s)-2,3-dihydroxypropyl]aoxy}(hydroxy)phosphoryl]oxy}-2-[(6e)-hexadec-6-enoyloxy]propyl (8e)-octadec-8-enoate was the single ligand associated with CTNNB1 that was used for virtual screening against the FDA and TCM libraries. Unfortunately, neither library yielded a small molecule candidate with greater than a 50% Tanimoto score. The ligand of PLCG2 did yield several interesting drug leads in both the FDA-approved library and the TCM library. Fludarabine has been tested in oral cancer cell lines and found to be effective in inducing cell apoptosis. 77 Ganciclovir an HIV drug has also been tested in oral cancer and found to have a clinical effect on cell differentiation. 78 Entecovir and didanosine are drugs used in the treatment of HepB infection and HIV and have not yet been tested on oral cancer. Small molecules in the TCM database meeting shape overlay filters were [4-[2,6-dimethylmorpholin-4-yl)sulfonylphenyl]-[4-(2-phenoxyethyl)piperazin-1-yl]-methanone which has not yet been used on oral cancer cell lines. There was overlap between drug leads for ligands of UBB and PLCG2. This is due to the structural similarity between cytosine and guanine ligands used as reference molecules for similarity matching. Drug leads for UBB included Cytarabine (Cancer), Fludarabine (Cancer), Azacitidine (myelodysplastic syndrome), Gemcitabine (Cancer), and Lamivudine (HIV). Gemcitabine is unique in that it is the only drug of those listed, which has been approved by the FDA for use in oral cancer patients. For HDAC3, there were no matches exceeding a Tanimoto score threshold of 50% of the reference molecule when using a library of FDA-approved drugs. The TCM database did yield a match with an extract from Mallotus philippensis, a member of the Euphorbiacae plant family. An extract of this plant known as Rottlerin has been found to inhibit growth of colon cancer and breast cancer cells. 79,80 Quercetin and Diosmetin were other phenols found in citrus that were also identified as matches meeting Tanimoto thresholds. There are no reports of the effect of quercetin or diosmetin in oral cancer. It should be noted that Quercetin and Rottlerin have been noted in literature as promiscuous ligands that are often found in natural product in silico screenings. 81 A recent study reported quercetin as the number one natural product in terms of number of occurrences within the database. 82 The aforementioned studies do show that Rottlerin is related to the metastatic potential and viability of colorectal cancer cells. 79 Caution and validation of results with existing literature or carefully designed follow-up experiments should always be pursued to justify the results of promising in silico analyses. These analyses enhance knowledge of genes influencing treatment response in oral cancer. This study describes how pathway, network, and chemical informatics analysis can be paired with literature review to identify drug leads for oral cancer treatment. Reference ligands associated with influential genes in the Royal Chemistry Society Protein Database are listed along with drug leads and their corresponding Tanimoto similarity scores in Table 4.

Discussion
A machine learning approach known as random forest was used to identify genes influencing oral cancer treatment response specific to the platinum-based chemotherapy treatment type and the non-platinum-based chemotherapy treatment type. This article emphasizes the benefits of integrating the results of this line of analyses with pathway, network, and chemical informatics analysis to identify promising gene targets and drug leads. Biological plausibility of these findings was highlighted with a review of existing literature supporting the findings for pathways, genes, and small molecules that our reported approach identified as influential in oral cancer. The results of this work identify pathways influencing treatment response in platinum-based chemotherapy users, non-users, and those common to both users and non-users. Network analysis via Cytoscape allowed for the identification of those influential genes within each treatment modality group within the context of interconnected gene networks. The utility of random forest was underscored in that in addition to pathways it also provides a rank to each gene in its influence on treatment response. This approach is a low-cost method of prioritizing gene targets and drug leads. These methods are validated in that the genes identified have been shown to be associated with cancer progression in oral cancers and other cancers. Several drug leads identified were also shown to be effective in inhibiting oral cancer cells and were reported to be in different phases of the drug approval pipeline. A possible criticism of the method outlined in this study is that there is uncertainty in the degree of trust that should be extended to random forest measures of gene influence and the inference of importance to the pathways in which "influential" genes reside. To further such criticisms, a point could be made that the Tanimoto threshold of >50% similarity could be perceived as low and the 50% difference in the molecules compared may prevent activity and may also be shown to have toxicity for a given disease state. Given such uncertainties, it may seem that the evidence supporting these methods is tenuous. This study recognizes these criticisms; however, the counterpoints must be made that gene influence is not observed in a single sample of the data, but rather in more than 20 000 permuted samples of the data in which the top-ranked genes were found to be more influential than thousands of other genes. The computational intensity provided in this study (20 000 trees and 320 tries at each node of each decision tree) provides justification of the trust provided for each gene influence value. In respect to the results yielded by chemical informatics analysis, it is important to note that the outlined chemical informatics method was able to identify gemcitabine as a drug that has been approved for use in oral cancer by the FDA. This method also identified Fludarabine and Ganciclovir which have both been reported as providing significant reduction in oral cancer cell line progression and viability.
Drug leads were identified in both FDA and TCM libraries, and the benefit in expanding the number of libraries is that it increases the probability of finding a match meeting the Tanimoto threshold of >50%. The negative aspect of adding libraries is that if computational resources are not planned for accordingly, then the amount of time required to screen against each reference molecule will scale upward with library size. The tools used in this study were all open source and freely available, and a limitation to the adoption of this pipeline is that tools and their dependencies are distributed across different R repositories that may or may not be kept up to date. Combining these tools into a single package that allows for the identification of both gene targets and drug leads may enhance the pace of drug discovery pipelines. We have shown in this study that random forest is well suited to datasets with small observations and high number of features. Gene targets that have been shown (through literature review) to be associated with treatment response and cancer progression were identified through this study's use of random forest analysis. Stratifying this analysis by the type of chemotherapy received allows for interpretation of influential genes and pathways within the context of treatment. Indeed, the lack of overlap in the importance of genes from one treatment modality to another highlights the idea that the gene expression patterns influencing platinumbased treatment response differ from those gene expression patterns influencing non-platinum treatment response.
It is likely that there is bias inherent to the stratification of patients by chemotherapy treatment type. Chemotherapy treatment is associated with clinical variables like clinical stage, tumor size, and tumor grade, as well as gender and socioeconomic quintile. [83][84][85] This study attempted to address these confounders by including them within the bag of randomly selected features available for construction of the random forest model. By integrating chemical informatics analyses, random forest results can be translated into lists of drug leads for each target gene. This method identified drug leads that have already entered or passed phase 3 trials. Our review of identified drug leads and comparison with existing annotations show that the chemical informatics methods described can identify small molecules with therapeutic potential. This study provides the impetus for further exploration of the role of the identified small molecules in oral cancer treatment response and the targeting of those genes identified as most influential by our series of analyses. This study also serves as a model for researchers identifying gene targets in rare cancers where the number of cases is limited.

Author Contributions
JS was the sole contributor to manuscript composition, experimental design, analysis, and production of all tables and figures in this manuscript.

Supplemental Material
Supplemental material for this article is available online.