Systems Genetics of Liver Fibrosis: Identification of Fibrogenic and Expression Quantitative Trait Loci in the BXD Murine Reference Population

The progression of liver fibrosis in response to chronic injury varies considerably among individual patients. The underlying genetics is highly complex due to large numbers of potential genes, environmental factors and cell types involved. Here, we provide the first toxicogenomic analysis of liver fibrosis induced by carbon tetrachloride in the murine ‘genetic reference panel’ of recombinant inbred BXD lines. Our aim was to define the core of risk genes and gene interaction networks that control fibrosis progression. Liver fibrosis phenotypes and gene expression profiles were determined in 35 BXD lines. Quantitative trait locus (QTL) analysis identified seven genomic loci influencing fibrosis phenotypes (pQTLs) with genome-wide significance on chromosomes 4, 5, 7, 12, and 17. Stepwise refinement was based on expression QTL mapping with stringent selection criteria, reducing the number of 1,351 candidate genes located in the pQTLs to a final list of 11 cis-regulated genes. Our findings demonstrate that the BXD reference population represents a powerful experimental resource for shortlisting the genes within a regulatory network that determine the liver's vulnerability to chronic injury.


Introduction
Liver fibrosis is a common consequence of chronic injury.Inducing agents vary from hepatotoxins, metabolic disorders and autoimmune reactions to viral infections.A characteristic feature of the fibrotic response is the ongoing repair mechanism resulting in an excessive accumulation of extracellular matrix [1,2].Fibrosis may progress to liver cirrhosis, which is characterized by severe distortion of liver architecture and impaired function.Of note, in patients with an exposure to similar environmental risk factors, the progression of liver fibrosis varies markedly.Based on the rate of fibrosis progression, patients may be classified as 'slow' or 'rapid fibrosers' [3].These interindividual differences have been attributed to environmental, but also to genetic (and epigenetic) factors [1,4,5].Several fibrogenic gene variants have been identified, e.g.complement component 5 (Hc) [6], other chemoattractants such as the chemokine CXCL9 [7] and the chemokine receptor CXCR3 [8] or metabolic enzymes like the triglyceride hydrolase adiponutrin (PNPLA3) [9,10].Moreover, two recent genome-wide association studies [11,12] identified a set of novel potential susceptibility genes for liver fibrosis, including PNPLA3, but no specific networks underlying fibrogenesis were reported.However, due to the large number of factors involved, the systematic identification of genetic determinants and networks affecting hepatic fibrosis remains a major challenge.
Systems genetics is a powerful method to dissect the underlying mechanisms of complex traits, including predisposing gene networks and environmental variants [13,14].The key experimental set-up is to make use of a genetic reference population.Here, we availed of the BXD set of recombinant inbred (RI) mouse lines, which are inbred progeny of F2 intercrosses of the inbred mouse strains C57BL/6J and DBA/2J [13,15].RI lines are especially suited as a mapping panel, since they form an immortalized set of isogenic lines [16,17], and a large number of animals and phenotypes per genome can be analyzed under standardized experimental conditions, thus lowering environmental noise.This improves the yield of information for the detection of genetic loci linked to trait variation, known as quantitative trait locus (QTL) mapping [13,15].In previous studies we have demonstrated that the parental strains C57BL/6J and DBA/2J show significant phenotypic variation of key fibrogenic parameters and therefore differ in their fibrosis susceptibility [6,7,18].Since these strains also vary in four million genetic sites across their genome [19], they provide the phenotypic and genetic diversity necessary for mapping studies in liver fibrosis.Furthermore, with more than 13,000 genetic markers and over 3,000 phenotypic records the BXD lines are one of the best-characterized murine reference panels [13,20,21].
Our aim was to determine new gene variants that affect hepatic fibrosis and to apply a systems genetics approach for the identification of gene networks that are critical for fibrosis phenotypes.Therefore, we characterized differences in fibrosis susceptibility of BXD lines after induction of liver fibrosis with carbon tetrachloride (CCl 4 ) and generated toxicogenomic, hepatic expression profiles by microarray analyses.Afterwards, we associated the genetic variation in our population with transcript variation in order to identify determinants of gene expression in liver fibrosis.In addition to single QTL and gene-gene interaction studies, the combination with expression genetics provided novel insights into potential networks modifying hepatic fibrogenesis.

Animals and experimental design
C57BL/6J, DBA/2J, B6D2 F1 hybrids and BXD lines were obtained from The Jackson Laboratory (Bar Harbor, ME) or from Oak Ridge Laboratory (lines BXD43, BXD51, BXD61, BXD62, BXD65, BXD68, BXD69, BXD73, BXD75, BXD87, BXD90), and were bred in the facility of the Neurobsik consortium from the VU University Amsterdam.The mice were maintained in a mouse facility under controlled environmental conditions.
In addition to the parental strains and the F1 hybrids, we studied 35 BXD lines with an average of six mice per sex and line, resulting in a total of 581 mice.Liver fibrosis was induced at eight weeks of age.For a period of six weeks, CCl 4 was administered by intraperitoneal injections twice weekly (0.7 mg CCl 4 /kg body weight in mineral oil, final volume 50 ml).Forty-eight hours after the last injection, the animals were anaesthetized with isofluran and killed by cervical dislocation.Instantly, blood was collected from the vena cava inferior, and tissue samples of liver and spleen were harvested.Liver samples were divided into five separate lobes.Whole liver weight, spleen and body weights were noted.
The animal studies were conducted according to all relevant welfare regulations and the Animal Care and Use Committee for Saarland University approved the protocols (TV Nr. 10/2008).

Phenotypic characterization of hepatic fibrosis
We measured the following quantitative CCl 4 -induced phenotypes: hepatic collagen contents (hydroxyproline levels and collagen areas) as quantitative measures and fibrosis stage as semiquantitative measure in histological liver sections.In 35 BXD lines, C57BL/6J, DBA/2J and F1 hybrids, collagen contents were determined in liver hydrolysates from snap frozen specimens of the right hepatic lobe.The assay is based on photometric measurement of the collagen specific amino acid hydroxyproline (Hyp) and follows the slightly revised protocol of Jamall et al. [18,22].
For the histological assessment of liver injury, formalin-fixed left lobes (4%, v/v) were available from 29 BXD lines, strains C57BL/ 6J and DBA/2J as well as B6D2 F1 hybrids.Each lobe was cut into 3-4 cross sections and embedded in one paraffin block.To detect collagen fibers, paraffin sections were stained with Sirius red [18].The staging of fibrosis was performed using a semiquantitative scoring system adapted from the system of Batts and Ludwig [18,23,24], principally differentiating the stages F0 to F4 ('F-scores').
Furthermore, stained collagen areas were quantified by morphometric analysis, using a semiautomatic system for image analysis (Stingray F146C IRF Medical camera, K0 type progressive scan CCD, Germany, and HistoQuant image morphometry software, 3DHistech, Budapest, Ungary).Mean collagenous areas (mm 2 ) were calculated by setting a threshold capturing Sirius red stained areas of collagen.One representative field (magnification 100x) was chosen from each liver section (avoiding arteries of a diameter .100mm), and the mean percentage of the stained area to whole area (field of vision) was calculated.Liver injury was assessed by serum alanine aminotransferase (ALT) activity.After CCl 4 challenge, blood was collected in a terminal procedure as described above.Blood was centrifuged for 20 min with 20006g at 4uC.Serum was diluted with 0.9% (v/v) NaCl, and ALT levels were determined in the central laboratory of Saarland University Medical Center according to the IFCC reference method (Cobas, Roche Hitachi, Indianapolis, IN) [25].

Microarray analysis of hepatic expression profiles
Total RNA was isolated from snap frozen individual liver samples (,30 mg) of 30 BXD strains, the parental strains and B6D2 F1 hybrids, using the RNeasy mini kit (Qiagen, Hilden, Germany).Three female mice per strain were analyzed after CCl 4 treatment for six weeks as described above, resulting in a total of 99 liver samples.
RNA quality was verified by measurement of the RNA integrity number (2100 Bioanalyzer, Agilent, Santa Clara, CA).Whole genome profiles of the fibrotic livers were performed using Gene Chip Mouse Gene 1.0 ST arrays (Affymetrix, Santa Clara, CA).For the normalization of robust multi-array average (RMA) intensity estimates of each transcript, RMAs were transformed into log 2 -values.Then the data of each single array was converted to their Z-scores, so that each array has a mean of 0 and a standard deviation of 1. Subsequently all values were multiplied by 2, and a value of 8 was added.Accordingly, all final values are positive and all datasets have an average of 8 units.One unit of expression on this scale resembles approximately a two-fold difference in expression.
The Affymetrix expression data set comprised a total of 34,760 records, which were assigned to 22,349 annotated genes.We determined the mean gene expression values for each probe set.As internal control, the strain distribution patterns of eQTLs with a Mendelian (monogenic) expression pattern (e.g.Alad, deltaaminolevulinate dehydratase; Hc, hemolytic complement; Tceanc2, transcription elongation factor A) were determined to show a perfect match to those of their closest markers, verifying that there were no errors of strain assignment in this data set (not shown).
All normalized transcript data are available in the GeneNetwork database (accession ID GN325, database name SUH BXD Liver CCl4-treated Affy Mouse Gene 1.0 ST (Jun11) RMA).
GeneNetwork is an open-access database that collates genomic information of diverse experimental crosses and reference panels as well as phenotypic data from miscellaneous research groups [26].

Statistics
Data generation, statistical analysis and graph creation were performed with SPSS Statistics 21 (IBM, Ehningen, Germany).As appropriate, mean and median values were further used for QTL analysis.Phenotypic robustness for each strain was assessed by the standard errors of the means.Mean and median values were trimmed by identifying and omitting outliers after graphical inspection of the data in box plots.Trait values of each BXD line were analyzed for each sex separately, as well as for the combined data sets of female and male mice.All fibrosis trait data were uploaded into the GeneNetwork database (accession IDs 14355-14396).Pearson's correlation was used to correlate fibrosis data among themselves and to BXD phenotypes.
Heritability (h 2 ) of the fibrosis traits was calculated as the ratios of between-strain variances to within-strain variances in sexspecific data sets [27].Within-strain and between-strain variances were calculated with analysis of variance (ANOVA) [28].
pQTL analysis.The identification and mapping of phenotypic QTLs (pQTLs) was performed by linking trait data to genotypes at known genetic marker loci [13,29].All phenotypic data were integrated into the GeneNetwork database.For the identification of single QTLs, interval mapping analyses were performed across all chromosomes [30].The parental strains were included in all mapping analysis.For composite interval mapping (CIM), a single genetic marker representing an identified QTL region was included as covariate, increasing the power to identify QTLs on other chromosomes by removing the effect of the preeminent QTL [31].CIM was performed for every phenotype, choosing genetic markers with highest likelihood ratio statistics (LRS) at each single QTL region or interacting loci from pairwise interaction scans.
The significance of a hypothetical QTL was estimated from the LRS [30].Genome-wide significance was evaluated by testing 2,000 permutations [32], which specified a significant threshold corresponding to a genome-wide p-value (p G ) of 0.05, and a suggestive threshold corresponding to p G = 0.63.Confidence intervals of chromosomal regions spanning QTL positions were specified as 1.5 logarithm of the odds (LOD) support intervals.
eQTL analysis.In this study, we aimed to infer causal mechanisms for phenotypic variation.Observing effects on gene expression that result from variants in the identified genomic region increases the confidence that this locus harbors causal candidates underlying the phenotype.Therefore we followed a pre-defined selection strategy for candidate genes, herein we chose to restrict our candidate search to the genes located in significant and phenotype-overlapping pQTL regions (Figure S1 in File S1).Genes were identified using the GeneNetwork/UCSC Genome browser [33].Since gene expression is, at least in part, heritable, differences in mRNA expression levels over the panel of BXD lines can be used to map regulatory expression quantitative trait loci (eQTLs).Therefore hepatic expression levels were used as quantitative traits and implemented into interval mapping analyses to identify regulatory loci.An expression quantitative trait locus (QTL) for a specific transcript was denoted cisQTL if the associated marker was localized within a 10 Mb distance of the gene position [34].The respective gene was called a cis-regulated quantitative trait gene (cisQTG).
Integrated search for candidate genes.The lists of cisQTGs for the fibrosis phenotypes were refined following an explorative in silico data analysis (Figure S2 in File S1).We applied the following selection criteria, of which at least one had to be fulfilled by the cisQTG to qualify as a potential fibrogenic susceptibility gene: (i) significant Pearson's correlation coefficient (p,0.05corresponding to r.0.36 or r,20.36) between cisQTG and any of the fibrosis phenotypes (collagen area, Hyp, F-score); (ii) non-synonymous single nucleotide polymorphism (nsSNP) that differs between the parental strains B6 and D2.nsSNPs leading to amino acid substitutions in the coding regions were identified using the GeneNetwork variant browser; and (iii) differential hepatic regulation of cisQTGs in expression data sets after CCl 4 challenge and saline treated control livers.QTL regions were matched using QTLminer [35] and the GeneNetwork dataset for saline treated livers in females (accession ID GN312, database name GenEx BXD Sal Liver Affy M430 2.0 (Feb11) RMA Females).
Gene network construction.Pairwise correlation estimates of the expression values of designated candidate genes were calculated as Pearson's correlation coefficient and presented in a circular network graph using Cytoscape software [36].For purposes of illustration the intra-chromosomal correlations of genes on the same chromosome were omitted, and only interchromosomal correlations are shown in the network.Edge colors represent positive or negative correlations with r.0.36 or r,20.36, and node sizes indicate the connectivity of the genes.
All candidate genes including the three fibrosis phenotypes were additionally illustrated in a QTL heatmap.This heatmap visualizes the p-values of regulatory loci identified by genomewide linkage analysis of the traits, computed on the basis of permuation tests (n = 1,000).

Liver fibrosis in parental strains
A significant (p,0.0001)accumulation of collagen was observed in all livers after six weeks of CCl 4 challenge in comparison to untreated mice.The mortality rate of mice after CCl 4 challenge was 8%.The parental inbred strain DBA/2J is more susceptible to liver fibrosis than the C57BL/6J strain, indicated herein by significantly (p,0.05)increased collagen areas and higher hepatic collagen (Hyp) contents; this is paralleled by liver injury as assessed by serum ALT activities (Figure 1).

Liver fibrosis phenotypes in the BXD reference panel
We noted marked differences in liver fibrosis among the BXD lines.Hepatic collagen contents varied widely (mean 6 SD 386.96141.5 mg Hyp/g liver; Figure 2A), and histopathological fibrosis scores correlated significantly (p,0.001) with hepatic collagen levels (Figure 2B).Furthermore, we observed significant (p,0.0001)line differences with respect to hepatic collagen levels and fibrosis scores (Figures 2C-E) as well as clinical-chemical parameters (ALT) of liver damage (Figure 2F).Semiquantitative fibrosis scores ranged from F1 (perivenular fibrosis) to F4 (pronounced fibrosis) (Figure 2G-K).
Heritability of liver fibrosis as determined by mean h 2 values was similar for all fibrosis traits, ranging from 0.5160.18(Hyp) to 0.5760.02(F-score) and 0.5960.01(collagen area).Overall, hepatic collagen area was identified as the most heritable trait (h 2 = 0.59).h 2 for hepatic Hyp concentrations showed a difference between male (0.36) and female mice (0.87), whereas h 2 for the other traits did not differ between sexes.

Genome-wide mapping of liver fibrosis phenotypes (pQTLs)
Single QTL genome scans identified 28 trait associated loci (Table S1 in File S2) with LRS scores above the suggestive threshold (see Methods) affecting liver fibrosis phenotypes.Among these significantly linked loci were detected by composite interval mapping (p G ,0.05) on chromosomes 4, 5, 7, 12, and 17 (Table 1, Figure S3 in File S1).For all loci except the QTL on chromosome 12, alleles of the fibrosis-susceptible strain DBA/2J increased the trait values (Table S1 in File S2).Five QTLs on chromosomes 2, 5, 7, 13 and 15 conferred susceptibility to more than one phenotype (Table S2 in File S2), whereas the loci on chromosomes 4, 12 and 17 were specific for a single phenotype.We found that 16 of 28 QTLs (46%) were sex-specific, i.e. they were only found in data from male or female mice (Table S1 in File S2).The remaining 11 loci were detected in the combined datasets; but showed significant effects either in male or female, consistent with a predominant phenotypic effect of one sex.In addition, the chromosome 2 QTL at 174.5-181.5 Mb was detected in all datasets tested.QTLs for hepatic Hyp levels were mainly based on female datasets; this was in line with the higher h 2 scores for this trait in female mice.For collagen area and F-score QTLs, no predominance of a single sex was observed.

Genome-wide mapping of fibrosis-associated eQTLs
Nine pQTL regions on chromosomes 2, 4, 5, 7, 12, 13, 15 and 17 were further dissected using eQTL mapping (Figure 3).The nine selected loci either conferred significant linkage or were associated to more than one phenotype.eQTLs were cisregulatory loci (cisQTLs) or trans-regulatory loci (transQTLs) (see Methods).By mapping cis-regulatory eQTLs within the nine most significant pQTL regions, we identified fibrosis-associated expression patterns that were locally regulated within the pQTLs.Table 1 summarizes the results of the eQTL analysis: On average, the pQTLs spanned an interval of 15.2 Mb, and these pQTL regions contained a total number of 1,351 annotated genes.The highest LRS score was observed for the QTL on chromosome 12, the only locus for which fibrosis susceptibility was conferred by the C57BL/ 6J allele.The QTL on chromosome 7 was the largest QTL with a high gene density.Overall, we identified 68 regulatory markers (eQTLs) within the pQTL regions.Thirty regulatory markers were identified as cisQTLs.Using the upper limit of suggestive thresholds for genome-wide significance as determined by permutation tests (LRS$12.0),these markers were linked to 85 genes in close proximity (,10 Mb).The associated genes are potentially cis-regulated genes (cisQTGs) in the pQTL regions.Of note, all cisQTLs also demonstrate trans-regulation of additional genes outside the pQTL regions (not shown).Applying this combined analysis of pQTLs and eQTLs, we reduced the number of potential fibrosis candidate genes from 1,351 to 85 (Table 2, Figure S2 in File S1).

Selecting fibrogenic candidate genes
In further selection steps, we inferred key regulatory candidates of fibrosis among the 85 cisQTGs.This strategy refined the list of 85 cisQTGs to 55 potential profibrogenic candidate genes that fulfilled at least one selection criterion (Table S3 in File S2).First, Pearson's correlation of cisQTG expression with any fibrosis phenotype (collagen area, Hyp, F-score) identified 30 significantly (p,0.05)correlated genes (indicated by dark gray boxes).Next, we determined covariance of hepatic mRNA expression patterns in unchallenged and CCl 4 challenged livers of the BXD reference lines, using QTLminer analysis.This analysis indicated that eQTLs could be distinguished into (A) fibrosis-specific cisQTLs that showed differential regulation between the basal state and after fibrosis induction, and (B) fibrosis-independent cisQTLs, i.e. the genes were cis-regulated in both groups.We speculate that differential regulation of class A genes in fibrotic livers identifies more relevant modifiers of fibrogenesis.In total, 45 cisQTGs were differentially regulated (class A), while 40 cisQTGs were cisregulated in both normal and fibrotic livers (class B).Finally, we identified 169 genes in the pQTL regions with nsSNPs that segregated between the two parental lines C57BL/6J and DBA/ 2J.Among the 85 cisQTGs, 22 have nsSNPs in coding regions.

Generating a fibrosis gene network
We generated a circular network graph (Figure 4), using the pairwise correlation estimates of hepatic expression levels for the 55 fibrosis candidate genes (Table S3 in File S2).The nodes represent the genes, which were arranged according to their chromosomal localizations, and edges show significant (p,0.05)inter-chromosomal correlations between gene expression levels.This unique network contains a total of 115 inter-chromosomal correlations; for reasons of simplicity the graph does not show the 93 intra-chromosomal correlations.The genes with the highest   The heatmaps represent significant interval mapping results on the indicated mouse chromosomes, separately for male and female mice as well as the combined data set (without inclusion of covariates); the QTL plots below illustrate composite interval mapping results (with 'background' QTLs as covariates, restricted to significant QTLs or overlapping loci for different phenotypes).Color coding of the heatmaps is as follows: Grey/black regions indicate the absence of genotype to phenotype linkage.Blue to green regions correspond to suggestive and significant linkage, respectively, with C57BL/6J alleles being associated with higher trait values.Red to yellow regions correspond to suggestive and significant linkage, respectively, with an association of DBA/2J alleles with higher values.doi:10.1371/journal.pone.0089279.g003 In addition, we illustrated the regulatory eQTLs of the candidate genes and the loci for fibrosis phenotypes (pQTL) in a QTL heatmap (Figure S4 in File S1).Regulatory gene clusters on chromosomes 2, 4, 5, 7 and 12 co-localized with loci of fibrosis phenotypes, as indicated by arrows in the QTL heatmap.The large eQTL on chromosome 7 showed two differentially regulated gene clusters: The genes located on distal chromosome 7 (63.1-74.2Mb) displayed additional regulatory loci on chromosome 5, which were absent for the genes between 51.2 and 53.5 Mb.

Discussion
Previously we identified individual loci that confer genetic susceptibility to hepatic fibrosis in different experimental crosses of inbred mice [6,7,18].Here we report the first systems genetics analysis of fibrosis in the BXD murine reference panel that allows the integration of multiple traits [13].The genetically mosaic BXD inbred lines display significant variation of quantitative fibrosis phenotypes, consistent with polygenic inheritance of liver fibrosis [1,4,5].By correlating phenotypes and known BXD genotypes in a genome-wide QTL analysis, we identified multiple pQTLs, nine with genome-wide significance and several with sex-dependent effects.Sex-specific differences are observed in various liver diseases such as (non-)alcoholic fatty liver diseases and hemochromatosis and might be due to sex hormone-regulated mechanisms or sex-specific gene variants.In addition to the phenotypic characterization of fibrosis, we generated a comprehensive expression dataset that represents the first genome-wide transcriptome analysis of hepatic fibrosis in a murine reference panel.
By stepwise bioinformatic analyses [37], we were able to reduce the number of 1,351 genes located in the nine pQTL regions to a set of 55 creedal fibrogenic candidate genes (Figure S2 in File S1).For this analysis, we focused our search on the genes that are cisregulated during fibrogenesis.To minimize the false discovery rate, the genes underwent a subsequent careful explorative analysis.They were considered to be relevant for fibrosis when their expression levels correlated with the fibrosis traits or showed differential regulation in healthy and fibrotic livers.In addition, we screened these cisQTGs for nsSNPs in the parental strains of the BXD panel, which might structurally and/or functionally affect protein functions.After all these steps, eleven genes fulfilled all selection criteria (Afm, Fan1, Hsd17b14, Klk1b5, Klk1b21, Klk1b22, Klk1b26, Napsa, Nomo, Nin, Susd1).For the majority of these genes there is no established connection to hepatic fibrosis or little information about their function in liver, although kallikreins exert known functions in the activation of inflammation, wound healing, and liver regeneration [38].In particular, we consider Afm as interesting candidate with a so far unknown role in hepatic fibrogenesis.Afm is a member of the albumin gene family that was shown to be differentially regulated by hepatocyte nuclear factors 1a and 1b in mice [39].Interestingly it functions as carrier of vitamin E [40], which has recently been reported to ameliorate liver fibrosis in fatty liver disease in mice and humans [41,42].A study by Kim et al. [43] revealed that afamin acts as a chemokine activating the Akt-signalling cascade, at least in osteoblasts.Because expression profiling showed Afm to be differentially expressed in mouse liver, this observation suggests similar regulatory effects across organs.
A further proof of principle is that several of the significant QTL regions include potential candidates that have previously been associated with fibrosis progression, in particular the chemokine ligand Cxcl10, Nr1h2 (nuclear receptor subfamily 1, group H, member 2, a.k.a.LXR), and Tnc (tenascin C).Cxcl10 encodes a chemokine that promotes hepatic inflammation by leukocyte Table 2. Summary of fibrosis-associated pQTL and eQTL regions.recruitment [44].Additional in vitro experiments in primary mouse hepatocytes detected a time-dependent induction of Cxcl10 expression levels after treatment with Tgfb1, supporting its role as profibrogenic candidate gene (R.H., R.L. and F.L.; unpublished observations).The nuclear receptor LXR might exert antifibrotic effects, since it was shown to reduce hepatic stellate cell activation, and therefore inhibit the production of profibrogenic cytokines [45].Tnc is an extracellular matrix glycoprotein expressed by hepatic stellate cells and myofibroblasts in liver, where it increases cytokine expression and ameliorates leukocyte transmigration [46,47].Whereas Tnc expression is absent in naı ¨ve livers, it is strongly induced during enhanced cell turnover as seen in wound repair [48,49].It also represents an endogenous ligand of Tlr4, which promotes innate immune responses during fibrogenesis [50,51].Further correlation network analysis of our transcriptome datasets showed that Tnc is highly interconnected with other fibrosis-associated genes (R.H. and F.L., unpublished observations); in particular, it is significantly correlated with the hepatic expression levels of collagens (Col1a2: r = 0.76; p = 0.002; Col3a1: r = 0.75; p = 0.001) and Tgfb1 (r = 0.74; p = 0.002).
To illustrate the complexity of fibrosis susceptibility, we combined the genes and their expression correlations in a large 'fibrosis network', indicating the inter-chromosomal connectivity of fibrogenic genes.In addition, the genetic architecture of coregulated gene networks were visualized by eQTL heatmaps, which show that the regulatory gene clusters of fibrosis susceptibility genes co-localize with pQTLs (Figure S4 in File S1).These findings reflect the selection concept for the candidate genes, with gene clusters representing networks with a causal relationship to fibrosis phenotypes.Gene repression studies (e.g.siRNA) might help to further dissect the directionality of the gene effects and fibrosis gene networks.Recently genome-wide association studies (GWAS) in patients with chronic liver diseases have identified profibrogenic gene variants [5].Genes associated with fibrosis progression in patients with chronic HCV infection were MERTK (c-mer proto-oncogene tyrosine kinase), RNF7 (ring finger protein 7), and TULP1 (tubby like protein 1) [12].In contrast, the genes PNPLA3 (adiponutrin), GCKR (glucokinase regulator) and TRIB1 (tribbles homolog 1) were associated with fibrosis phenotypes in non-alcoholic fatty liver disease [52], with PNPLA3 demonstrating the most consistent effects across studies [53].Although our expression dataset shows that all murine orthologs are expressed in fibrotic mouse liver (mean expression scores 7.2-10.4;R.H. and F.L., unpublished observations), solely Pnpla3 was located in a pQTL on chromosome 15.However, it was apparently not cis-regulated in our panel and did not display nsSNPs in the parental strains.The lack of cross-species validation might be due to the distinct induction of fibrosis in our model, the specific charcteristics of the patients included in the GWAS (viral hepatitis, fatty liver disease), different phenotypic effect sizes across species, or non-polymorphic regions within the BXD panel.
This study illustrates how advances in the methodologies of systems genetics with the use of a murine reference panel lead to the identification of a potential disease network for liver fibrosis.The BXD lines represent an appropriate reference population with phenotypic segregation of fibrosis phenotypes due to genetic variation.Our findings indicate that it is essential not to focus on single fibrogenic QTLs, but on gene clusters as modifiers of fibrosis susceptibility.Although several existing confounders are being controlled for within this approach, further developments in gene mapping and functional validation will contribute to the translation of our experimental findings to patients with liver fibrosis.

Figure 1 .
Figure 1.Phenotypic characterization of parental strains after six weeks of CCl 4 injections.Liver fibrosis was assessed by morphometric (A) and biochemical (B) measurement of hepatic collagen (Hyp) contents.Hepatic inflammation was measured by serum ALT activities (C).Sirius red staining of hepatic collagen showed circumferential fibrosis in C57BL/6J mice (D) and pronounced fibrosis in DBA/2J mice (E), corresponding to mean F-scores of 2.060.1 and 3.960.1,respectively.doi:10.1371/journal.pone.0089279.g001

Figure 3 .
Figure3.QTLs for hepatic fibrosis in the BXD murine reference population.The heatmaps represent significant interval mapping results on the indicated mouse chromosomes, separately for male and female mice as well as the combined data set (without inclusion of covariates); the QTL plots below illustrate composite interval mapping results (with 'background' QTLs as covariates, restricted to significant QTLs or overlapping loci for different phenotypes).Color coding of the heatmaps is as follows: Grey/black regions indicate the absence of genotype to phenotype linkage.Blue to green regions correspond to suggestive and significant linkage, respectively, with C57BL/6J alleles being associated with higher trait values.Red to yellow regions correspond to suggestive and significant linkage, respectively, with an association of DBA/2J alleles with higher values.doi:10.1371/journal.pone.0089279.g003

Figure 4 .
Figure 4. Fibrosis network graph generated by correlating the hepatic expression levels of the candidate genes.The graph presents the inter-chromosomal correlations of 51 genes, except for the four candidate genes Adamts17, Gm9860, Klk22 and Ogfr with exclusively intrachromosomal correlations.Node color and shape illustrate the chromosomal localization of the gene.The size of each node indicates the degree of connectivity, with larger nodes having higher number of correlated genes.The edges show Pearson's correlation coefficients (r) as follows: solid green lines: r.0.5; solid red lines: r,20.5;dotted green lines r.0.36; dotted red line r,20.36.doi:10.1371/journal.pone.0089279.g004

Table 1 .
Chromosomal regions of pQTLs with significant genome-wide LRS values determined by single QTL scans and CIM.