Epidemic dispersion of HIV and HCV in a population of co-infected Romanian injecting drug users

Co-infections with HIV and HCV are very frequent among people who inject drugs (PWID). However, very few studies comparatively reconstructed the transmission patterns of both viruses in the same population. We have recruited 117 co-infected PWID during a recent HIV outbreak in Romania. Phylogenetic analyses were performed on HIV and HCV sequences in order to characterize and compare transmission dynamics of the two viruses. Three large HIV clusters (2 subtype F1 and one CRF14_BG) and thirteen smaller HCV transmission networks (genotypes 1a, 1b, 3a, 4a and 4d) were identified. Eighty (65%) patients were both in HIV and HCV transmission chains and 70 of those shared the same HIV and HCV cluster with at least one other patient. Molecular clock analysis indicated that all identified HIV clusters originated around 2006, while the origin of the different HCV clusters ranged between 1980 (genotype 1b) and 2011 (genotypes 3a and 4d). HCV infection preceded HIV infection in 80.3% of cases. Coincidental transmission of HIV and HCV was estimated to be rather low (19.65%) and associated with an outbreak among PWID during detention in the same penitentiary. This study has reconstructed and compared the dispersion of these two viruses in a PWID population.


Introduction
The use of contaminated needles is an effective way of spreading HIV; it accounts for 10% of infections worldwide and for more than 40% in some regions and countries (Eastern Europe, South America, and East and Southeast Asia) [1,2]. Seroprevalence of HCV antibodies is very PLOS ONE | https://doi.org/10.1371/journal.pone.0185866 October 9, 2017 1 / 17 a1111111111 a1111111111 a1111111111 a1111111111 a1111111111 frequent among people who inject drugs (PWID) rising up to 90%; in the general HIV infected population it raises only to about 10% [3][4][5]. HIV-HCV co-infection affects the natural history of both diseases: increased mortality risk, increased AIDS and liver-related disease mortality, accelerated hepatic fibrosis that leads to cirrhosis and/or hepatocellular carcinoma [6]. Both viruses have a high level of genetic diversity within and between hosts and a high rate of evolution. For HIV-1, several groups, subtypes and circulating recombinant forms (CRFs) are described; for HCV seven different genotypes have been described and each genotype has distinct subtypes [7,8]. The high evolutionary rate allows us to study the transmission history of both viruses by analysing their genomic sequences with phylogenetic methods.
In recent years, an important HIV outbreak was reported among PWID in Romania, with 30.6% of the total newly diagnosed HIV cases in the country in 2012 [9]. The most affected region was Bucharest and its suburbs where more than 90% of Romanian PWID are living [10]. The Romanian PWID infected with HIV are mainly young males, unemployed, often diagnosed in penitentiary, who inject both heroin and new psychoactive drugs ("legal highs"), and are frequently co-infected with HCV. Two predominant HIV strains were identified in PWID: subtype F1 which had already been circulating in sexual and nosocomially infected patients and the newly endemic CRF14_BG, most probably originating in Spain or Portugal [11,12].
Previous studies established the role of phylogenetic analysis in assessing and characterizing transmission networks within particular populations, but these were mainly performed in specific populations of HIV patients. However, few studies have been performed that reconstruct simultaneously the transmission history of HIV and HCV viruses in the same population of HIV-HCV co-infected patients [13,14] and none analysed the overlap between HIV and HCV transmission networks. Our objective was to reconstruct the history and the transmission patterns of HIV and HCV, to assess the likelihood of simultaneous transmission of these two viruses in a population of co-infected PWID and to analyse the congruence between HIV and HCV transmission networks.

Study population
We have included in this study a number of 117 injecting drug users who were diagnosed with HIV infection between 2011 and 2014 in the National Institute for Infectious Diseases 'Matei Bals', the reference center for HIV genotyping. All patients were HCV co-infected and naive to antiretroviral treatment. A number of 35 non-PWID, HCV mono-infected patients, were also included in the study.
Blood samples and clinical and epidemiological data were collected through a questionnaire survey. Laboratory testing (i.e. HIV viral load, CD4 count, HCV serology, HCV viral load) was performed with in vitro diagnostic (IVD) tests.

Ethics statement
The study was conducted according to the Declaration of Helsinki and approved by the Ethical Committee of the National Institute for Infectious Diseases 'Matei Bals'. All the subjects included in this study were adults who provided their written informed consent.

HIV genotyping
Total viral RNA was extracted from 1 ml of plasma using the automated NucliSens Easy-MAG nucleic acid extraction system (BioMerieux, Marcy-l'Étoile, France), according to the manufacturer's instructions and eluted in 25μl of water. Ten microliters of RNA were used for HIV genotyping. Reverse transcription and amplification of the HIV-1 pol gene (1301 nucleotides, HXB2: 2253-3554) were performed using the Viroseq HIV-1 Genotyping System (Celera Diagnostics, Alameda, CA) and sequenced with the ABI 3500 Genetic Analyzer (Applied Biosystems, Foster City, CA). The primary data were analysed using Sequencing Analysis Software Version 3.7 (Life Technologies) and the generated sequences were assembled with ViroSeq 2.8 HIV-1 Genotyping System Software (Celera Diagnostics, Alameda, CA). REGA HIV-1&2 automated subtyping tool version 3.0 was used to assign the HIV-1 subtype [15]. Unclassified samples were further analysed in order to identify possible recombination breakpoints using Simplot v3.5.1 software and the Los Alamos HIV-1 subtype reference dataset (sliding window: 400-nt, T:t ratio = 2.0, model of evolution: Kimura twoparameter, bootstrap: 1000 replicates). All sequences were screened for hypermutation using the Hypermut 2.0 algorithm [16].
The NS5b and NS3 amplicons were bidirectionally sequenced using BigDye1 Terminator system v 1.1 (Life Technologies) and ABI 3500 Genetic Analyzer (Applied Biosystems). Seqscape version 2.7 (Applied Biosystems) was used to assemble and generate consensus sequences. Genotyping was done using the publicly available algorithm Oxford HCV Automated Subtyping Tool (Version 2.0) [19].

Phylogenetic, phylodynamic analysis and divergence times of transmission networks
Reference sequences were selected from public databases [20, 21] taking into account genotype / subtype, genomic region and geographic origin. To retrieve control sequences for transmission cluster reconstruction, the most similar HIV sequences were retrieved from our and Los Alamos database using BLAST. Duplicates and clones were discarded. Of 421 selected control sequences, those representative to geographical and time dispersion of the patients sample, were subsequently used for each analysis. Controls for all the HIV-1 subtypes and CRFs identified among Romanian PWID were selected as follows: CRF14_BG (Spain, Portugal, Romania), subtype F1 (Angola, Romania), subtype B (European countries), CRF35_AD (Afghanistan, Iran). CRF14BG epidemics from Spain-Portugal [11] and Greece [12] were previously linked to the Romanian HIV epidemic among PWID. Subtype F1 strains from Romanian PWID were previously shown [11] to originate from strains already circulating in Romania and formerly described as related with strains from Angola [22,23].
A similar procedure was used for the HCV dataset, where a number of 304 control sequences were used. Multiple alignments were generated using Muscle [24].When available, concatenate alignments of NS5b and NS3 sequences were used in the Bayesian analyses.
Phylogenetic analyses were performed using both maximum-likelihood (ML) and Bayesian approaches, after selecting the best-fit nucleotide substitution model using JModelTest (GTR +4Γ+I for HIV-1 and TN93 for HCV). Datasets for each HIV-1 subtype and HCV genotype identified among Romanian PWID were constructed and analyzed separately. ML trees were inferred using FastTree software [25] and Bayesian phylogenetic analyses were performed with BEAST version 1.8.1 [26] using a Bayesian Skyline coalescent tree prior [27] and the uncorrelated lognormal relaxed clock model [28]. Two Markov Chain Monte Carlo (MCMC) runs were computed separately for 10 7 generations with a burn-in of 10%. The output of the MCMC analysis was tested for convergence by means of effective sampling size (ESS>100) using the program Tracer v1.6 [29]. The maximum clade credibility (MCC) tree was summarized with TreeAnnotator program [30].
Transmission clusters were identified as clusters in the ML phylogenetic tree with SH-like support higher or equal to 0.9 and were further confirmed with Bayesian phylogenetic analyses if supported by a Bayesian posterior probability value of 1.0 [31,32]. The temporal signal of both datasets was evaluated with TempEst [33]. Data were statistically evaluated and plotted using GraphPad Prism v6.0. The Mann-Whitney test was performed to identify differences between groups and the Chi-square or Fisher exact test were used to evaluate the distribution in particular subgroups; p-values <0.05 were considered significant.

Characteristics of the population
The population of 117 co-infected PWID included in this study were diagnosed with HIV between 2011 and 2014. The diagnosis of HCV and HIV was simultaneous for 71 patients (60.7%), while for the remaining 46 individuals HCV diagnosis was established earlier than the HIV diagnosis. Most of the PWID (64.1%, 75) were asymptomatic at the HIV diagnosis, while 35.9% were diagnosed in more advanced HIV infection stages, classified as stage B (27%, 32) and C (9%, 10) according to CDC and WHO Staging Systems [34]. The 35 control HCV mono-infected patients were diagnosed during the same period of time, also in the National Institute for Infectious Diseases 'Matei Bals'.

Disease severity differs when comparing patients infected with HIV-1 subtype F1 and with CRF14_BG-like strains
HIV-1 subtype analysis showed that 67.5% (n = 79) of the patients were infected with subtype F1 viruses; recombinant forms present in a significant proportion (26%; n = 30) were CRF14_BG (n = 19, 16.24%) and unique recombinant forms (URF) of CRF_14BG and F1 (n = 11, 9.4%). The PWID infected with these recombinant forms had significantly lower CD4 counts at baseline than patients infected with subtype F1 strains. It is true, though, that an estimation of the time to the infecting event was not possible in either case. No significant differences between these two groups were observed when HIV viral load (VL) or CDC stages were analysed (Fig 1). However, a higher percentage of patients infected with CRF14_BG were diagnosed at stages B and C than the group infected with F1 (47% vs. 34%). Few patients (n = 8, 6.8%) were infected with other types of HIV: subtype B, CRF35_AD and URFs between subtype A1 or B and F1 (Table 1).

HCV viral loads are higher in patients with HIV co-infection
All 117 patients were tested to assess the HCV genotype by analysing NS5b and/or NS3 regions of HCV. Both NS5b and NS3 sequences were available for 110 patients, while for the remaining patients only the NS5b region was available. The genotype assessment indicated that the majority of PWID were infected with genotype 1a (n = 47, 40%) and 1b (n = 40, 34%), but genotypes 3a and 4 (a, d) were also identified ( Table 1). HCV subtype assignment was concordant in the two genomic regions in the majority of cases (n = 108, 98.2%), with only two exceptions of genotype 1 sequences (2749bh2013 and 2727bh2013) assigned as 1a in NS3 and 1b in NS5b. Unlike HIV-HCV co-infected PWID, the individuals in the HCV mono-infected group were exclusively infected with genotype 1b strains. Significant differences (Fig 1) in HCV VL  were observed when comparing mono-infected with co-infected patients, with higher viral loads in the group of co-infected.

All identified HIV clusters are large and originated around 10 years ago
Phylogenetic analysis of Romanian PWID using control sequences from Romania and other parts of the world showed that 90 HIV sequences (76.9%) were part of well-supported phylogenetic clusters, suggesting high levels of transmission networking. Three major HIV local transmission networks were identified by phylogenetic analysis, two within subtype F1 and one with CRF14_BG (Fig 2). One of the subtype F1 clusters included 29 PWID patients, 24 men As opposed to HIV where the virus has spread in only three large clusters, HCV dispersion occurred through multiple transmission networks. In total, 97 HCV sequences (82.9%) clustered in 13 transmission clusters. In comparison, within the HCV mono-infected control group only 8 out of 33 (24.2%) HCV sequences were found to form phylogenetic clusters (χ 2 = 42.18, p<0.001).
Four well supported clusters of different sizes were observed in both genotypes 1a and 1b; the largest consisted of 23 sequences (genotype 1a) and 13 sequences (genotype 1b) ( Table 2). Three small genotype 3a clusters were observed, while genotype 4a and 4d sequences identified among PWID formed each one single well supported cluster (Fig 3).
Molecular clock analysis of HCV transmission clusters revealed important differences in origin and starting time points. The differences in divergence time were observed within as well as among genotypes. Based on relaxed clock and constant size population model implemented in a Bayesian MCMC framework, the oldest PWID network belongs to genotype 1b and originated around 1983 (95%HPD: 1968-1994), followed by the larger HCV genotype 1a cluster, dated back to 1987 (Figs 3 and 4). However, the time to the most recent common ancestor (tMRCA) estimated for a sub-group of this 1a cluster consisting in 22 out of 23 sequences was more recent (2004( , 95%HPD 1999.
Genotype 1b was found to be the first introduced in the studied population: the oldest cluster was estimated to be 30 years old in PWID and more than 50 years in non-PWID. Genotypes 3 and 4 are the most recently introduced in this risk group population, dating back 2-3 years before sampling. The divergence times for all important PWID networks, HIV and HCV, are presented in Table 2 and graphically illustrated in

Coincidental transmission of HIV and HCV is uncommon and occurred in 19.65% of the patients
We have comparatively analysed the patterns of patients distribution in HIV and HCV transmission clusters (Fig 5). Eighty patients were present both in HCV and in HIV transmission clusters (68.4% of all co-infected patients; 82.5% and 88.9% of all HCV and HIV patients present in clusters, respectively). Furthermore, 70 patients (59.8%) were present in the same HIV and HCV cluster with at least one other patient, indicating transmission of HIV and HCV through (at least partly) overlapping transmission chains.
Phylogenetic analysis and, in particular, measuring the difference in divergence times between viral lineages within a co-infected individual has been proven useful in assessing the frequency of viral coincidental transmission in PWID in the absence of seroconversion data [35]. We have used the calculated difference of mean divergence times (ΔtMRCA) between HIV and HCV lineages within a patient to estimate the time period lapsed between the acquisition of the two viruses. Low ΔtMRCA indicates a higher likelihood that the two viruses might have been co-transmitted to a susceptible individual, while higher values of ΔtMRCA suggest that the viruses were transmitted at different and distinct time points. For the majority of the patients this analysis suggested that HCV was introduced first, followed by HIV. However, for 5 patients, HIV was introduced first; all of them were infected with subtype F1 and genotype 3a or 4 HCV strains. Twenty three patients had ΔtMRCA between HIV and HCV lineages of less than one year, suggesting potential co-transmission of HIV and HCV. These 23 patients (19.65%) were present in the same HIV and HCV cluster with at least one other patient, suggesting a potential simultaneous transmission (co-transmission) of both infections. We also found that such transmission events are largely associated with particular HCV transmission clusters, specifically clusters c and l (Fig 5). Going back to the records of these patients, we found that twenty-two of those 23 patients were in detention at the moment of study enrolment. This population has a statistically significant higher probability of coincidental transmission (p<0.0001) when compared to the rest of the patients. Five of these patients were infected with CRF14_BG and the rest with subtype F1; the infective HCV genotypes were 1 or 4 (Fig 5).

Discussion
The PWID represent no more than an estimated 0.2-0.5% of the world's population, but they account for approximately 5-10% of all people living with HIV and 6.8% of the persons infected with HCV [36,37]. Since 2008 injected amphetamine type stimulants ("legal highs") became available on the market in Romania and were low-priced. As a result, the number of PWID increased dramatically in recent years so that Romania reported an HIV outbreak among PWID that started in 2011; more than 90% of PWID being concentrated in Bucharest and its surroundings [10,38]. Numerous people who were injecting heroin have switched to the new class of drugs due to the lower price and because they were legal at that time. An important behavioral characteristic is that the latter require up to 10 injections every day. This fact, which was coincident with an  interruption in needle exchange programs, greatly impacted the needle-sharing behavior in the country [39]. Preliminary analyses indicated the circulation of CRF14_BG strains in PWID, HIV strains associated with lower CD4 counts and X4 tropism at baseline, as compared with the subtype F1 infected PWID [11,40].
In this study we have identified, characterized and compared the viral transmission networks of HIV and HCV in a population of HIV-HCV co-infected PWID-using a group of HCV mono-infected patients as controls-and we evaluated the time of introduction of these viruses in this population and explored the probability of co-transmission of the two infections.
The clinical and epidemiological impact of HIV-HCV co-infections in PWID has been described by several studies [6,41], but only one analysed in parallel the evolutionary history of these two viruses within this particular risk population [35]. Viral genetic information expressed as nucleotide sequences and analysed with phylogenetic and population genetics has been shown to be very useful in identifying and characterizing transmission clusters of rapidly evolving viruses such as HIV and HCV [42,43].
Phylogenetic analysis performed in this study showed high levels of transmission networking among Romanian PWIDs: three HIV clusters and thirteen HCV clusters have been identified. The majority of HIV strains fell in two F1 clusters that were nested within Romanian reference strains. This suggests that F1 strains circulating for a long time among Romanian patients, infected sexually and nosocomially were the origin of the subtype F1 outbreaks in PWID and that these strains were restricted to Romania. The CRF14_BG strains circulating in PWID are related with those described earlier in Spain and Portugal in the same risk group population [11]. However, due to the limited number of CRF14_BG sequences available in public databases, direct filiation between the two epidemiologic events cannot be ascertained. Previous results indicated that the CRF14_BG sequences from Romanian PWID were closely related with the Greek PWID sequences and formed a single monophyletic cluster. Molecular clock analysis revealed that CRF14_BG strains circulating in PWID from Athens originated from one or several strains circulating in Bucharest [12]. Our phylodynamic analyses indicated that all HIV transmission networks originated at similar time points, dating back 7 to 10 years ago, coincident with the availability of these new injectable amphetamine type stimulants on the market. The patients infected with CRF14_BG had significantly lower CD4 counts at baseline than those infected with F1 subtype strains despite the fact that phylogenetic and seroconversion data suggest their being infected at similar time points. Although we have included in our analysis all the PWID presented in our clinical settings during this period of time, the number of patients infected with CRF14_BG recombinants was smaller than those infected with subtype F1 and this might have impact on the statistical analysis. However, these data are in agreement with previous studies showing an association between CRF14_BG infection, CD4 T cells decline and disease severity, characterized by CXCR4 usage [11,40,44].
On the other hand, regardless of the HIV infecting strain, HCV viral loads were significantly higher in co-infected PWID than in HCV mono-infected patients. This difference does not seem to be related to the HCV genotype. Moreover, when comparing HCV mono-infected patients with co-infected PWID carrying HCV genotype 1b, the difference in viral load continued to be statistically significant (S1 Fig). Several studies have shown an accelerated disease progression in HIV-HCV co-infected patients, HIV infection increasing the HCV replication, hepatic inflammation and impairing HCV-specific immune response [6]. Since CRF14_BG induced a massive depletion of CD4 T cells it was expected to be associated with impaired HCV immune control and thus higher replication and hepatic fibrosis. However, no significant difference in HCV VL was observed between CRF14_BG and subtype F1 infected patients or between HCV VL and CD4 count. Possibly, longer infection periods could be needed to impair the established HCV-specific immune response, since in most of the cases, HCV infection preceded HIV infection.
Previous studies have shown the genotype 1b to be mainly acquired through blood transfusion, while genotypes 1a and 3a are associated with intravenous drug use [45]. In other European countries, mostly HCV genotype 4 has been identified in PWID: subtype 4a in Greece, 4d in Italy and The Netherlands and 4c and 4d in Spain were linked to drug usage [46][47][48][49].
Our results indicate the circulation of multiple HCV genotypes in co-infected Romanian PWID: 1a, 1b, 3a, 4a and 4d. However, only HCV genotype 1b was found in the HCV mono- infected group, indicating a higher complexity of the HCV epidemic among co-infected PWID.
According to phylodynamic analyses, HCV genotype 1b was the first to be introduced in Romania, followed by genotype 1a. The oldest 1b cluster was estimated to be 30 years old in PWID and more than 50 years in non-PWID. Previous studies indicated that the circulation of this genotype in Romania was almost exclusive for a long period of time [50,51]. Genotypes 3 and 4 are much more recent and linked to intravenous drug usage; with an estimated origin around 10 years ago, they might also be linked to the availability of injected amphetamine type stimulants [52].
The differences in transmission clustering between HCV and HIV might be explained in part by their estimated acquisition time. Previous observations based on classical epidemiological data indicated that HCV was introduced in PWID much earlier, being associated with heroin use while HIV infection is linked with new psychoactive drugs [11]. Phylodynamic analysis presented in this study sustains this observation, showing that HIV strains were introduced in PWID in recent years whereas HCV genotypes 1a and 1b had been circulating several decades before. However, HCV transmission networks of genotypes 3 and 4 apparently started later, as they were dated as very recent (2010-2011), thus indicating a potential linkage of these genotypes with new psychoactive drugs.
When comparing the transmission chains of HIV and HCV, we evaluated both the congruence between HIV and HCV transmission chains, as well as the likelihood of simultaneous transmission of both infections. We found that most patients (68.4%) were present both in HIV and HCV transmission clusters and that for the majority of the patients' transmission of HIV and HCV occurred through (at least partly) overlapping transmission chains (Fig 5). We used the method reported by Ng TK and collaborators [35] to evaluate the frequency of viral coincidental transmission. Applying this method to Romanian co-infected PWID, we found that, for most of the patients, and in accordance with epidemiological data, HCV infection was acquired before HIV.
Nevertheless, we identified 23 cases of possible co-transmission (19.65%), corresponding to patients with ΔtMRCA<1 and clustering in partly overlapping transmission chains. Interestingly, these transmission events were associated with particular clusters of genotype 1a and 4a. Further investigation lead us to the finding that 22 of these 23 patients were in detention in the same penitentiary. While largely overlapping transmission chains would be expected, given the fact that this population forcely shares the same contact network while in detention, a higher likelihood of simultaneous transmission of the infections is not necessarily expected in this context. Therefore, our findings call for further investigation about what type of behavioral patterns explain the high rate of co-transmission of HIV and HCV in this penitentiary when compared to the rest of the PWID population. One potential bias in our analysis could be the subevaluation of possible super-infections and co-transmissions: phylogenetic analysis based on population sequencing considers only the predominant viral infective strain. Coincidental transmissions could be underestimated in our analysis if the transmitted strains did not become visible within the viral population.
The HIV epidemic among Romanian PWID continues to represent a major problem and HCV co-infection contributes to the worse prognostic within this specific risk population. This study clarifies and compares the patterns of spread and transmission of these two viruses in Romanian PWID, suggesting that a higher rate of transmission of HIV rather than HCV is associated with behavioral changes mainly following the availability of new types of injected drugs ('legal highs'). While overlapping of HIV and HCV transmission chains is common, the estimated transmission of both viruses within a short time span is less common, with a rate of only 19.65% in this population, and largely associated with drug injection within the same penitentiary. Our results illustrate the need for better prevention policies in the PWID population, in particular and more urgently in the prison population, where an outbreak of simultaneous HIV and HCV transmissions seems to be occurring.
One limitation of the study is that the inmate population is overrepresented in the PWID so that extrapolation to the general population should be cautious. It should be though noted that, within this rather small group, no significant clinical or epidemiologic differences were noted according to the inprison vs. outprison status.

Accession numbers
86 HIV nucleotide sequences reported in this paper were newly deposited in GenBank and have the following accession numbers: KX158994-KX159079. The remaining 31 HIV sequences were previously published in Genbank [11] and have the following accession num-