- Email: [email protected]
Response to comment on “Toxicogenomics in human health risk assessment”
Toxicology and Applied Pharmacology 236 (2009) 257–260
Contents lists available at ScienceDirect
Toxicology and Applied Pharmacology j o u r n a l h o m e p a g e : w w w. e l s ev i e r. c o m / l o c a t e / y t a a p
Response to Letter to the Editor Response to comment on “Toxicogenomics in human health risk assessment”
To the Editor: We thank Drs. Eric Carlson and Jay Silkworth for their letter regarding the two genomics papers we published in this journal recently [Toxicol Appl Pharmacol 231, 165–178 (2008) and Toxicol. Appl. Pharmacol 231, 179–196 (2008)]. We have carefully considered their concerns and surmise that they did not fully understand the purpose and content of the papers. Before we address the concerns specifically, we add a general comment that the focus of these two papers was not, as implied by Drs. Carlson and Silkworth, to use genomic data as a stand-alone endpoint in risk assessment and to point out that there was no mention of risk assessment in either paper. Our goal was to determine mode(s) of action for chemical-induced toxicities using genomics as a tool. This preliminary information will guide future studies using proteomics and be used to compliment existing neurochemical data. Genomics is a valuable tool for defining upstream events to morphological and physiological outcomes. Although we recognize the promise of genomics in human health risk assessment, we do not believe that it can be used as a stand-alone approach. Therefore, our aim is to use genomics in concert with proteomics, biochemical events, and structural and functional changes to support human health risk assessments that are scientifically sound. We have clearly outlined this approach as a schematic in the second paper (Scheme 1; Royland and Kodavanti, 2008). The five specific criticisms of Drs. Carlson and Silkworth are addressed below in the order in which they were raised. Their first concern was that the use of a single high dose of Aroclor 1254 in our study is not relevant to human health risk assessment. Although we used a single 6 mg/kg in this preliminary genomics study, several other studies conducted by our lab included a lower dose of 1 mg/kg, at which significant changes in several neurochemical endpoints were observed (Yang et al., 2003). In addition, as the Aroclor 1254 was administered to the dams, the exposure in utero and via lactation would result in a much lower dose to the pups. Further, we saw no decrease in weight gain in treated dams, no significant difference in pup mortality or in pup postnatal weight gain. The paper cited by Drs. Carlson and Silkworth described a 2-year cancer study in adult rats (Mayes et al., 1998) where the rats were dosed with PCBs (25 to 200 ppm) for 2 years. Our experiments, in contrast, are developmental studies where the pups were exposed to the compounds both in utero and through lactation from dosed dams during times that are known to be critical for the development of different brain regions. As for the relevance of the dose, we have previously determined the PCB levels in blood, brain and other tissues of pups from Aroclor 1254-treated dams and found levels comparable to those seen in occupationally exposed humans (WHO, 1993). PCBs are ubiquitous chemicals that bio-accumulate in the environment. Therefore, humans potentially can be exposed from a number of 0041-008X/$ – see front matter. Published by Elsevier Inc. doi:10.1016/j.taap.2009.01.020
sources such as drinking water and food products (fish, meat, and dairy products). In selected regions PCB levels have been reported to range from 100 to 450 ng/l in drinking water and at levels of about 200 mg/kg fresh weight in food products (WHO, 1993). Further, Bush and Kadlec (1995) reported that zebra mussels from the Niagara River showed PCB levels in the mg/kg range. In humans, PCB levels in the adipose tissue of occupationally exposed workers have been reported to range from 2.2 to 290 ppm (WHO, 1993) which encompasses the range of 289 ppm detected in the adipose tissue of PND 21 pups perinatally dosed with 6 mg/kg Aroclor 1254 in the present study (Yang et al., 2003). Likewise, blood concentrations of PCBs detected in capacitor manufacturing workers were as high as 3.5 μg/ml (Wolff, 1985), which is similar to the level of 1.44 μg/ml found in the PND-21 pups in the present study. Thus, the doses used in our present genomic studies and in our earlier neurochemical, morphological, and functional measures would result in tissue levels in rats consistent with those found at the high end of the general population or in occupational workers. The second concern raised was regarding the relevance of the developmental neurotoxicity studies of PCBs. There are numerous studies from scientists all over the world reporting PCB-induced neurotoxicity in multiple animal models as well as in epidemiological studies in humans (see review, Kodavanti, 2005). A PubMed search with the key words “PCBs and neurotoxicity” identified more than 130 publications, including several review articles. While it is true that the neurobehavioral effects in rats reported by some investigators were subtle and no dose-response was observed, neurobehavioral data from mice and monkeys consistently showed significant differences (Schantz et al., 1989; Eriksson and Fredriksson, 1996). In addition, neurophysiological, morphological, and neurochemical data from rats as well as mice reproducibly showed PCB-induced effects (see review, Kodavanti 2005). It is clear from this extensive literature that PCBs interfere with neurodevelopment when exposure occurs during the rapid brain growth spurt. However, it is not clear if the effects on neurodevelopment are the only critical effects for PCB exposure. Drs. Carlson and Silkworth also raised relevance issues based on the differences in the brain developmental windows between rodents and humans. Because of this difference in developmental windows between rodents and humans, our dosing regimen covered the entire perinatal developmental window (in utero and lactational periods) to better relate rodent studies to the effects seen in humans. It is well known that the rat pup developmental state at birth corresponds approximately to the third trimester in humans; consequently lactational exposure in rat would be equivalent to exposure during the last third of in utero development in humans. Therefore the information obtained from our rodent studies, which cover the entire developmental window, will be useful to understand the effects of these chemicals on humans. Regarding the third concern, we are fully aware of the advantages of multiple testing correction to reduce the number of false positives and did test the data both ways. However, it is not uncommon to publish without including this relatively stringent step (for example, see Wei et al., 2008 and Lee et al., 2008). In cases where high
258
Response to Letter to the Editor
stringency is not required (i.e., in analysis where one is more interested in all possible players rather than trying to isolate the biggest player or in research vs. clinical studies), where expression differences tend to be small (as is the rule in the brain with over 200 cell types) or where small sample size (in our case 3 arrays/condition) can result in few genes being identified as significantly different, one is justified in using the less stringent analysis. To improve the analysis, our data were optimized for low expression genes and to minimize the number of possible false positives. The number of analyzed genes was reduced by removing any transcript for which the raw expression value was not more than 20 above background in at least one condition. This resulted in approximately one third fewer total transcripts analyzed, removed potentially questionable data points, and reduced significantly the statistically determined number of potential false positives. The data were then corrected for poor PM/ MM probe pairs with the Reduction of Invariant Probes (REDI) analysis before normalizing using GC Robust Multiarray Analysis (GC RMA), which is reported to improve discrimination of low expression genes (Wu et al. 2003). Filtering on a fold change after identifying differentially expressed genes by statistical analysis (i.e. t-test) further increases the stringency of the analysis and makes it more likely that the identified genes are truly differentially expressed and not found by chance alone. As to the 1.5 fold cutoff, the analysis algorithms for genomic data (e.g., GC RMA Normalization and REDI Probe Level Enhancement) have greatly improved, allowing for more accurate and lower discrimination of significant differences. A two fold change is no longer considered necessary for expression difference significance. This is particularly true in the brain where there is a vast variety of cell types (even within explicit brain areas) and expression differences are often quite low. In the brain, fold changes less than 2.0 are often cited (e.g. Dong, et al., 2005; Hakak, et al., 2001; Kaiser et al., 2004). Statistical significance was determined before calculating fold change, therefore, all listed genes were different at p ≤ 0.05. Also in regard to the third concern, we disagree that the RT-PCR data support the occurrence of false positives. In the developmental study (Royland et al., 2008, Figure 9), 10 of 12 comparisons in the cerebellum between array and RT PCR match quite closely. The two remaining apparent mismatches had, in one case (Cacnb1, control), no statistical difference between PND7 and PND14 for either protocol and for the second case (Gabrb3, Aroclor 1254), up-regulation for both techniques, but more with the RT-PCR. The results are less optimal for the hippocampus, but in each case where there is an apparent mismatch, the RT-PCR value is greater — not an unusual occurrence for the more sensitive PCR technique. That the PCR picked up some differences not found on the arrays refutes rather than supports the idea of false positives with the array analysis. Likewise, the RT-PCR data for the mode of action paper (Figure 10) closely matches the array data in most cases. One probe set against Cdh10 gave different values between the two methods at PND7, but matched well at PND14. Similarly the Gabrb3 data is an apparent mismatch in the hippocampus, but not in the cerebellum. We do not have a definitive explanation for the apparent discrepancies, but, as stated in the paper, they could be due to differential isoform expression at the different ages and/or tissues. Drs. Carlson's and Silkworth's claim that too much was made of the significance of the Aroclor 1254-induced changes in Gabrb3 is not warranted. Gabrb3 was mentioned only in the Materials and methods and/or Results sections of the manuscripts. In the Materials and methods, discussion was limited to the reasoning behind selecting it for RT-PCR confirmation (Royland et al., 2008). In the Results sections it was limited to two places: 1) to the results of the RT PCR studies (Royland et al., 2008; Royland and Kodavanti, 2008), where it was one of 2 transcripts identified in the array study to be differentially expressed in Aroclor 1254-treated animals in both cerebellum and hippocampus at PND7 (Royland and Kodavanti, 2008) and 2) to when citing literature for all the RT PCR genes examined (Royland and Kodavanti, 2008), we included the statement
that “Literature reports suggest that GABA receptor activation may have roles in neural cell proliferation, migration and differentiation (Behar et al., 1996; Behar et al., 2000; Behar et al., 1999; Fiszman et al., 1999; Mejia-Gervacio and Marty, 2006)”. By preference, discussion in both papers centered on groupings of functionally-related genes in order to generate hypotheses for future studies and individual genes were only discussed in regards to correlating this study with literature reports. As to the fourth concern, we also take exception to the comment that we “fail[ed] to properly report and interpret the output of their Ingenuity™ pathway analysis”. While it's true we did not report the significance values for the networks displayed, this is a common omission. Network “scores” are not generally reported in the literature, perhaps because it would be meaningless to display a network that wasn't highly relevant from the multiple ones presented. Network scores for this analysis ranged from a high of 47 to a low of 12. These scores are based on a p-value calculation on the likelihood that network eligible molecules are part of the network by random chance and are equal to the negative exponent of this calculation. Thus our lowest score would indicate a p-value b10− 12. Unlike the p-values available for canonical functions in an alternate Ingenuity analysis function, p-values for network canonical overlays are not provided. It is also true as Drs. Carlson and Silkworth stated that “Ingenuity™ pathways are based on manually annotated databases derived from scientific literature”. The canonical pathways are over-laid on networks that start with our list of significantly differentially expressed genes (for example the 2 out of 14 genes mentioned for Aryl Hydrocarbon Receptor (AhR) signaling in Figures 5 and 6) but include genes from multiple studies in multiple tissues under multiple conditions and using different platforms — that is all possible connections reported in the literature. In addition, reported canonical pathway overlays were not chosen by selection (i.e. from predetermined functions of interest) or assigned randomly, but are those that demonstrated the most interactions among the numerous canonical functions available in the Ingenuity™ software, therefore most relevance, to the networks obtained. Thus it is not surprising that connections to genes are made that do not show up as significantly different in our study. The importance of the networks is to show possible interactions, report how our identified gene list relates to the literature and provide avenues for future investigations. The value of the networks presented is that they independently verified interactions with the AhR and calcium signaling, mechanisms reported in the literature to be important in PCB action, as well as those of less well reported functions in neuronal development (e.g. Axonal guidance) and synaptic function (e.g. synaptic LTP). Our interpretations of the Ingenuity™ data analyses were confirmed by the technical representatives at Ingenuity Systems™. In further regard to the fourth concern, the criticism that we incorrectly stated that changes in Ah receptor signaling were affected by Aroclor 1254 treatment is also incorrect. It is known that some components of the Aroclor 1254 such as mono-ortho PCBs to a lesser extent and non-ortho PCBs to a greater extent, affect AhR signaling (Safe, 1994). Non-ortho components made up a relatively small proportion of the total mix (0.02%) in our Aroclor 1254 (Kodavanti et al., 2001), therefore it is not surprising that their effects on AhRmediated expression were minimal. It was not clear as to what point Drs. Carlson and Silkworth are making in their comments regarding the genes, Nfia and Col1a2. The expression of these genes has been reported as being modulated by the AhR signaling pathway, as is that of Cyp1a2 gene the authors were interested in. In fact, nuclear factor I/ A (Nfia) is a transcription factor that binds, among other things, Cyp1a2, and is thus an integral part of that AhR-regulated pathway. In addition, Aroclor 1254 has been reported to modulate other nuclear transcription factors in the liver resulting in the up-regulation of the related gene, Cyp1a1 (Borlok et al., 2002). That we did not identify the Cyp1a2 gene in these studies may be due to the preponderance of non-
Response to Letter to the Editor
AhR-interacting congeners in the Aroclor 1254 mixture or due to the lack of Cyp expression in the brain. Very few studies on the presence of Cyp's in the brain have been reported and these generally report very low constitutive levels. Of the over 800 PubMed references on the Cyp1a2 gene, only 15 studies were in brain and only 3 of those investigated the AhR signaling pathway. The connection of the Col1a2 gene is via its regulation by TNF (Umannová et al., 2007), a well known AhR-regulated gene. Changes in expression of Col1a2, with its roles in cell migration, proliferation and adhesion, suggest an AhRmediated developmental effect. In regards to the ‘fifth’ criticism, it is indeed true that the same data set was analyzed differently in the two papers. However we fail to see how this is interpreted as a “flaw”. Genomic data sets are routinely analyzed and re-analyzed to pull out different types of information and these data are published separately. The difference here is that we chose to present the data in two companion manuscripts as the most efficient way to portray different aspects of this multi-dimensional study. The pathway analyses in the mode of action manuscript (Royland and Kodavanti, 2008) support and expand the developmental aspects outlined in the developmental manuscript (Royland et al., 2008). In the first manuscript (Royland et al., 2008) we evaluated global developmental effects by determining statistically significant changes in both constitutive gene expression and in gene expression following Aroclor 1254 treatment in two brain areas between PND7 and PND14. Genes in the individual condition lists are subsets of the total number of genes statistically different with age (see Figure 1, Venn diagram). In each case the top 100 differentially expressed genes were chosen as the most likely to be representative of the most impacted functional groups (complete gene lists are available from the journal on line as supplements). As it was our intent to examine functional correlates, we portrayed the relative number and direction of change for statistically different genes (PND7 to PND14) within functional groups for each condition and tissue (Figures 4 and 7). No statistics were necessary here as these lists are of genes identified that changed with age only in control or only in Aroclor-treated animals, but not in both. Likewise, the GO ontology analyses report p-values for each tissue/treatment and can be compared as to relative importance, but are not designed to report across treatment differences here. In any developmental study, developmental effects will be inseparable from treatment effects, thus it is not surprising that many genes were found in common between the two treatments. A similar functional analysis was carried out on these genes (Figures 3 and 6) to add greater depth to the study. Commentary on how the relative expression level was impacted for individual genes that were common between the conditions (Tables 2 and 4) was merely to highlight a possible general effect of Aroclor 1254 on gene expression. As stated in the figure legends and as is most often done in heatmap generation, the heatmaps display expression changes relative to the mean of all expression across all samples/conditions that are included in the heatmap (note, each figure contains heatmap representations of three individual clustering analyses). We believe all the information required to understand Figures 5 and 8 is included in the accompanying legends. We also believe that our analyses across developmental ages in this first manuscript (Royland et al., 2008), which is different from the more traditional comparison of control vs. treatment in the companion manuscript, is a novel approach that adds new insight into global Aroclor 1254 effects and its impact on neurodevelopment, providing possible new avenues for research. The concerns Drs. Carlson and Silkworth raised regarding control/ Aroclor comparisons were addressed in the companion manuscript on possible modes of action (Royland and Kodavanti, 2008). During the brain growth spurt occurring during much of the dosing regimen, individual brain areas are rapidly changing, each on its own timetable. This is apparent in the multivariate PCA analysis (Figure 9) showing that gene expression is variable by brain area, age and treatment (especially at PND7). These rapidly changing gene environments
259
resulted in an insufficient number of genes (see Venn diagrams, Figures 1 and 4) to do the kind of gene-by-gene pathway analysis across both time and dose as recommended by Drs. Carlson and Stillworth. However, we strongly believe that our grouping of transcripts into functional categories is appropriate. Individual gene expression is known to be affected by timing, homeostatic mechanisms, feedback loops, and a host of normal physiological processes. Reporting on the pathways that were enriched (i.e., with multiple transcripts involved) under our experimental conditions provides important information on possible modes of action. Despite the concerns raised by Drs. Carlson and Silkworth, we strongly believe that, although the current genomic data are preliminary, the results from this study are very useful in guiding complementary studies using proteomics and neurochemical measures to determine mode(s) of action for these chemicals. References Behar, T.N., Li, Y.X., Tran, H.T., Ma, W., Dunlap, V., Scott, C., Barker, J.L., 1996. GABA stimulates chemotaxis and chemokinesis of embryonic cortical neurons via calcium-dependent mechanisms. J. Neurosci. 16, 1808–1818. Behar, T.N., Scott, C.A., Greene, C.L., Wen, X., Smith, S.V., Maric, D., Liu, Q.Y., Colton, C.A., Barker, J.L., 1999. Glutamate acting at NMDA receptors stimulates embryonic cortical neuronal migration. J. Neurosci. 19, 4449–4461. Behar, T.N., Schaffner, A.E., Scott, C.A., Greene, C.L., Barker, J.L., 2000. GABA receptor antagonists modulate postmitotic cell migration in slice cultures of embryonic rat cortex. Cereb. Cortex 10, 899–909. Borlok, J., Dangers, M., Thum, T., 2002. Aroclor 1254 modulates gene expression of nuclear transcription factors: implications for albumin gene transcription and protein synthesis in rat hepatocyte cultures. Toxicol. Appl. Pharmacol. 181, 79–88. Bush, B., Kadlec, M.J., 1995. Dynamics of PCBs in the aquatic environment. Great Lakes Research Reviews 1, 24–30. Dong, H., Wade, M., Williams, A., Lee, A., Douglas, G.R., Yauk, C., 2005. Molecular insight into the effects of hypothyroidism on the developing cerebellum. Biochem. Biophys. Res. Commun. 330, 1182–1193. Eriksson, P., Fredriksson, A., 1996. Developmental neurotoxicity of four orthosubstituted polychlorinated biphenyls in the neonatal mouse. Environ. Toxicol. Pharmacol. 1, 155–165. Fiszman, M.L., Borodinsky, L.N., Neale, J.H., 1999. GABA induces proliferation of immature cerebellar granule cells grown in vitro. Brain Res. 115, 1–8. Hakak, Y., Walker, J.R., Li, C., Wong, W.H., Davis, K.L., Buxbaum, J.D., Haroutunian, V., Fienbert, A.A., 2001. Genome-wide expressin analysis reveals dysregulation of myelination-related genes in chronic schizophrenia. PNAS. 98, 4746–4751. Kaiser, S., Foltz, L.A., George, C.A., Kirkwood, S.C., Bemis, K.G., Lin, X., Gelbert, L.M., Nisenbaum, L.K., 2004. Phencyclidine-induced changes in rat cortical gene expression identified by microarray analysis: implications for schizophrenia. Neurobiol. Dis. 16, 220–235. Kodavanti, P.R.S., 2005. Neurotoxicity of persistent organic pollutants: possible mode(s) of action and further considerations. Dose Response 3, 273–305. Kodavanti, P.R.S., Kannan, N., Yamashita, N., Derr-Yellin, E.C., Ward, T.R., Burgin, D.E., Tilson, H.A., Birnbaum, L.S., 2001. Differential effects of two lots of Aroclor 1254: congener-specific analysis and neurochemical end points. Environ. Health Perspect. 109, 1153–1161. Lee, M.H., Kim, M., Lee, B.H., Kang, K.S., Kim, H.L., Yoon, B.I., Chung, H., Kong, G., Lee, M.O., 2008. Subchronic effects of valproic acid on gene expression profiles for lipid metabolism in mouse liver. Toxicol. Appl. Pharmacol. 226, 271–284. Mayes, B.A., McConnell, E.E., Neal, B.H., Brunner, M.J., Hamilton, S.B., Sullivan, T.M., Peters, A.C., Ryan, M.J., Toft, J.D., Singer, A.W., Brown, J.F., Jr., Menton, R.G., Moore, J.A., 1998. Comparative carcinogenicity in Sprague-Dawley rats of the polychlorinated biphenyl mixtures Aroclors 1016, 1242, 1254, and 1260. Toxicol. Sci. 41, 62–76. Mejia-Gervacio, S., Marty, A., 2006. Control of interneuron firing pattern by axonal autoreceptors in the juvenile rat cerebellum. J. Physiol. 571, 43–55. Royland, J.E., Kodavanti, P.R.S., 2008. Gene expression profiles following exposure to a developmental neurotoxicant, Aroclor 1254; pathway analysis for possible mode(s) of action. Toxicol. Appl. Pharmacol. 231, 179–196. Royland, J.E., Wu, J., Zawia, N.H., Kodavanti, P.R.S., 2008. Gene expression profiles in the cerebellum and hippocampus following exposure to a neurotoxicant, Aroclor 1254: developmental effects. Toxicol. Appl. Pharmacol. 231, 165–178. Safe, S., 1994. Polychlorinated biphenyls (PCBs): environmental impact, biochemical, and toxic responses, and implications for risk assessment. Crit. Rev. Toxicol. 24, 87–149. Schantz, S.L., Levin, E.D., Bowman, R.E., Heironimus, M.P., Laughlin, N.K., 1989. Effects of perinatal PCB exposure on discrimination-reversal learning in monkeys. Neurotoxicol. Teratol. 11, 243–250. Umannová, L., Zatloukalová, J., Machala, M., Krcmár, P., Májková, Z., Hennig, B., Kozubík, A., Vondrácek, J., 2007. Tumor necrosis factor-alpha modulates effects of aryl hydrocarbon receptor ligands on cell proliferation and expression of cytochrome P450 enzymes in rat liver “stem-like” cells. Toxciol. Sci. 99, 79–89. Wei, Y., Liu, Y., Wang, J., Tao, Y., Dai, J., 2008. Toxicogenomic analysis of the hepatic effects of perfluorooctanoic acid on rare minnows (Gobiocypris rarus). Toxicol. Appl. Pharmacol. 226, 285–297.
260
Response to Letter to the Editor
Wolff, M.S., 1985. Occupational exposure to polychlorinated biphenyls (PCBs). Environ. Health Perspect. 60, 133–138. World Health Organization, 1993. Environmental health criteria 140: polychlorinated biphenyls and terphenyls, In: Dobson, S., van Esch, G.J. (Eds.), nternational Programme on Chemical Safety, 2nd ed. World Health Organization, Geneva, Switzerland. Wu, Z., Irizarry, R.A., Gentleman, R., Murillo, F.M., Spencer, F., 2003. A Model Based Background Adjustment for Oligonucliotide Expression Arrays. Technical Report. John Hopkins University, Dept. Biostatistics Working Papers, Baltimore. Yang, J.-H., Derr-Yellin, E.C., Kodavanti, P.R.S., 2003. Alterations in brain protein kinase C isoforms following developmental exposure to a polychlorinated biphenyl mixture. Mol. Brain Res. 111, 123–135.
Joyce E. Royland Prasada Rao S. Kodavanti* U.S. Environmental Protection Agency, Neurotoxicology Division, B105-06, Research Triangle Park, NC 27711, USA E-mail address: [email protected] (P.R.S. Kodavanti). *Corresponding author. 27 January 2009
Contents lists available at ScienceDirect
Toxicology and Applied Pharmacology j o u r n a l h o m e p a g e : w w w. e l s ev i e r. c o m / l o c a t e / y t a a p
Response to Letter to the Editor Response to comment on “Toxicogenomics in human health risk assessment”
To the Editor: We thank Drs. Eric Carlson and Jay Silkworth for their letter regarding the two genomics papers we published in this journal recently [Toxicol Appl Pharmacol 231, 165–178 (2008) and Toxicol. Appl. Pharmacol 231, 179–196 (2008)]. We have carefully considered their concerns and surmise that they did not fully understand the purpose and content of the papers. Before we address the concerns specifically, we add a general comment that the focus of these two papers was not, as implied by Drs. Carlson and Silkworth, to use genomic data as a stand-alone endpoint in risk assessment and to point out that there was no mention of risk assessment in either paper. Our goal was to determine mode(s) of action for chemical-induced toxicities using genomics as a tool. This preliminary information will guide future studies using proteomics and be used to compliment existing neurochemical data. Genomics is a valuable tool for defining upstream events to morphological and physiological outcomes. Although we recognize the promise of genomics in human health risk assessment, we do not believe that it can be used as a stand-alone approach. Therefore, our aim is to use genomics in concert with proteomics, biochemical events, and structural and functional changes to support human health risk assessments that are scientifically sound. We have clearly outlined this approach as a schematic in the second paper (Scheme 1; Royland and Kodavanti, 2008). The five specific criticisms of Drs. Carlson and Silkworth are addressed below in the order in which they were raised. Their first concern was that the use of a single high dose of Aroclor 1254 in our study is not relevant to human health risk assessment. Although we used a single 6 mg/kg in this preliminary genomics study, several other studies conducted by our lab included a lower dose of 1 mg/kg, at which significant changes in several neurochemical endpoints were observed (Yang et al., 2003). In addition, as the Aroclor 1254 was administered to the dams, the exposure in utero and via lactation would result in a much lower dose to the pups. Further, we saw no decrease in weight gain in treated dams, no significant difference in pup mortality or in pup postnatal weight gain. The paper cited by Drs. Carlson and Silkworth described a 2-year cancer study in adult rats (Mayes et al., 1998) where the rats were dosed with PCBs (25 to 200 ppm) for 2 years. Our experiments, in contrast, are developmental studies where the pups were exposed to the compounds both in utero and through lactation from dosed dams during times that are known to be critical for the development of different brain regions. As for the relevance of the dose, we have previously determined the PCB levels in blood, brain and other tissues of pups from Aroclor 1254-treated dams and found levels comparable to those seen in occupationally exposed humans (WHO, 1993). PCBs are ubiquitous chemicals that bio-accumulate in the environment. Therefore, humans potentially can be exposed from a number of 0041-008X/$ – see front matter. Published by Elsevier Inc. doi:10.1016/j.taap.2009.01.020
sources such as drinking water and food products (fish, meat, and dairy products). In selected regions PCB levels have been reported to range from 100 to 450 ng/l in drinking water and at levels of about 200 mg/kg fresh weight in food products (WHO, 1993). Further, Bush and Kadlec (1995) reported that zebra mussels from the Niagara River showed PCB levels in the mg/kg range. In humans, PCB levels in the adipose tissue of occupationally exposed workers have been reported to range from 2.2 to 290 ppm (WHO, 1993) which encompasses the range of 289 ppm detected in the adipose tissue of PND 21 pups perinatally dosed with 6 mg/kg Aroclor 1254 in the present study (Yang et al., 2003). Likewise, blood concentrations of PCBs detected in capacitor manufacturing workers were as high as 3.5 μg/ml (Wolff, 1985), which is similar to the level of 1.44 μg/ml found in the PND-21 pups in the present study. Thus, the doses used in our present genomic studies and in our earlier neurochemical, morphological, and functional measures would result in tissue levels in rats consistent with those found at the high end of the general population or in occupational workers. The second concern raised was regarding the relevance of the developmental neurotoxicity studies of PCBs. There are numerous studies from scientists all over the world reporting PCB-induced neurotoxicity in multiple animal models as well as in epidemiological studies in humans (see review, Kodavanti, 2005). A PubMed search with the key words “PCBs and neurotoxicity” identified more than 130 publications, including several review articles. While it is true that the neurobehavioral effects in rats reported by some investigators were subtle and no dose-response was observed, neurobehavioral data from mice and monkeys consistently showed significant differences (Schantz et al., 1989; Eriksson and Fredriksson, 1996). In addition, neurophysiological, morphological, and neurochemical data from rats as well as mice reproducibly showed PCB-induced effects (see review, Kodavanti 2005). It is clear from this extensive literature that PCBs interfere with neurodevelopment when exposure occurs during the rapid brain growth spurt. However, it is not clear if the effects on neurodevelopment are the only critical effects for PCB exposure. Drs. Carlson and Silkworth also raised relevance issues based on the differences in the brain developmental windows between rodents and humans. Because of this difference in developmental windows between rodents and humans, our dosing regimen covered the entire perinatal developmental window (in utero and lactational periods) to better relate rodent studies to the effects seen in humans. It is well known that the rat pup developmental state at birth corresponds approximately to the third trimester in humans; consequently lactational exposure in rat would be equivalent to exposure during the last third of in utero development in humans. Therefore the information obtained from our rodent studies, which cover the entire developmental window, will be useful to understand the effects of these chemicals on humans. Regarding the third concern, we are fully aware of the advantages of multiple testing correction to reduce the number of false positives and did test the data both ways. However, it is not uncommon to publish without including this relatively stringent step (for example, see Wei et al., 2008 and Lee et al., 2008). In cases where high
258
Response to Letter to the Editor
stringency is not required (i.e., in analysis where one is more interested in all possible players rather than trying to isolate the biggest player or in research vs. clinical studies), where expression differences tend to be small (as is the rule in the brain with over 200 cell types) or where small sample size (in our case 3 arrays/condition) can result in few genes being identified as significantly different, one is justified in using the less stringent analysis. To improve the analysis, our data were optimized for low expression genes and to minimize the number of possible false positives. The number of analyzed genes was reduced by removing any transcript for which the raw expression value was not more than 20 above background in at least one condition. This resulted in approximately one third fewer total transcripts analyzed, removed potentially questionable data points, and reduced significantly the statistically determined number of potential false positives. The data were then corrected for poor PM/ MM probe pairs with the Reduction of Invariant Probes (REDI) analysis before normalizing using GC Robust Multiarray Analysis (GC RMA), which is reported to improve discrimination of low expression genes (Wu et al. 2003). Filtering on a fold change after identifying differentially expressed genes by statistical analysis (i.e. t-test) further increases the stringency of the analysis and makes it more likely that the identified genes are truly differentially expressed and not found by chance alone. As to the 1.5 fold cutoff, the analysis algorithms for genomic data (e.g., GC RMA Normalization and REDI Probe Level Enhancement) have greatly improved, allowing for more accurate and lower discrimination of significant differences. A two fold change is no longer considered necessary for expression difference significance. This is particularly true in the brain where there is a vast variety of cell types (even within explicit brain areas) and expression differences are often quite low. In the brain, fold changes less than 2.0 are often cited (e.g. Dong, et al., 2005; Hakak, et al., 2001; Kaiser et al., 2004). Statistical significance was determined before calculating fold change, therefore, all listed genes were different at p ≤ 0.05. Also in regard to the third concern, we disagree that the RT-PCR data support the occurrence of false positives. In the developmental study (Royland et al., 2008, Figure 9), 10 of 12 comparisons in the cerebellum between array and RT PCR match quite closely. The two remaining apparent mismatches had, in one case (Cacnb1, control), no statistical difference between PND7 and PND14 for either protocol and for the second case (Gabrb3, Aroclor 1254), up-regulation for both techniques, but more with the RT-PCR. The results are less optimal for the hippocampus, but in each case where there is an apparent mismatch, the RT-PCR value is greater — not an unusual occurrence for the more sensitive PCR technique. That the PCR picked up some differences not found on the arrays refutes rather than supports the idea of false positives with the array analysis. Likewise, the RT-PCR data for the mode of action paper (Figure 10) closely matches the array data in most cases. One probe set against Cdh10 gave different values between the two methods at PND7, but matched well at PND14. Similarly the Gabrb3 data is an apparent mismatch in the hippocampus, but not in the cerebellum. We do not have a definitive explanation for the apparent discrepancies, but, as stated in the paper, they could be due to differential isoform expression at the different ages and/or tissues. Drs. Carlson's and Silkworth's claim that too much was made of the significance of the Aroclor 1254-induced changes in Gabrb3 is not warranted. Gabrb3 was mentioned only in the Materials and methods and/or Results sections of the manuscripts. In the Materials and methods, discussion was limited to the reasoning behind selecting it for RT-PCR confirmation (Royland et al., 2008). In the Results sections it was limited to two places: 1) to the results of the RT PCR studies (Royland et al., 2008; Royland and Kodavanti, 2008), where it was one of 2 transcripts identified in the array study to be differentially expressed in Aroclor 1254-treated animals in both cerebellum and hippocampus at PND7 (Royland and Kodavanti, 2008) and 2) to when citing literature for all the RT PCR genes examined (Royland and Kodavanti, 2008), we included the statement
that “Literature reports suggest that GABA receptor activation may have roles in neural cell proliferation, migration and differentiation (Behar et al., 1996; Behar et al., 2000; Behar et al., 1999; Fiszman et al., 1999; Mejia-Gervacio and Marty, 2006)”. By preference, discussion in both papers centered on groupings of functionally-related genes in order to generate hypotheses for future studies and individual genes were only discussed in regards to correlating this study with literature reports. As to the fourth concern, we also take exception to the comment that we “fail[ed] to properly report and interpret the output of their Ingenuity™ pathway analysis”. While it's true we did not report the significance values for the networks displayed, this is a common omission. Network “scores” are not generally reported in the literature, perhaps because it would be meaningless to display a network that wasn't highly relevant from the multiple ones presented. Network scores for this analysis ranged from a high of 47 to a low of 12. These scores are based on a p-value calculation on the likelihood that network eligible molecules are part of the network by random chance and are equal to the negative exponent of this calculation. Thus our lowest score would indicate a p-value b10− 12. Unlike the p-values available for canonical functions in an alternate Ingenuity analysis function, p-values for network canonical overlays are not provided. It is also true as Drs. Carlson and Silkworth stated that “Ingenuity™ pathways are based on manually annotated databases derived from scientific literature”. The canonical pathways are over-laid on networks that start with our list of significantly differentially expressed genes (for example the 2 out of 14 genes mentioned for Aryl Hydrocarbon Receptor (AhR) signaling in Figures 5 and 6) but include genes from multiple studies in multiple tissues under multiple conditions and using different platforms — that is all possible connections reported in the literature. In addition, reported canonical pathway overlays were not chosen by selection (i.e. from predetermined functions of interest) or assigned randomly, but are those that demonstrated the most interactions among the numerous canonical functions available in the Ingenuity™ software, therefore most relevance, to the networks obtained. Thus it is not surprising that connections to genes are made that do not show up as significantly different in our study. The importance of the networks is to show possible interactions, report how our identified gene list relates to the literature and provide avenues for future investigations. The value of the networks presented is that they independently verified interactions with the AhR and calcium signaling, mechanisms reported in the literature to be important in PCB action, as well as those of less well reported functions in neuronal development (e.g. Axonal guidance) and synaptic function (e.g. synaptic LTP). Our interpretations of the Ingenuity™ data analyses were confirmed by the technical representatives at Ingenuity Systems™. In further regard to the fourth concern, the criticism that we incorrectly stated that changes in Ah receptor signaling were affected by Aroclor 1254 treatment is also incorrect. It is known that some components of the Aroclor 1254 such as mono-ortho PCBs to a lesser extent and non-ortho PCBs to a greater extent, affect AhR signaling (Safe, 1994). Non-ortho components made up a relatively small proportion of the total mix (0.02%) in our Aroclor 1254 (Kodavanti et al., 2001), therefore it is not surprising that their effects on AhRmediated expression were minimal. It was not clear as to what point Drs. Carlson and Silkworth are making in their comments regarding the genes, Nfia and Col1a2. The expression of these genes has been reported as being modulated by the AhR signaling pathway, as is that of Cyp1a2 gene the authors were interested in. In fact, nuclear factor I/ A (Nfia) is a transcription factor that binds, among other things, Cyp1a2, and is thus an integral part of that AhR-regulated pathway. In addition, Aroclor 1254 has been reported to modulate other nuclear transcription factors in the liver resulting in the up-regulation of the related gene, Cyp1a1 (Borlok et al., 2002). That we did not identify the Cyp1a2 gene in these studies may be due to the preponderance of non-
Response to Letter to the Editor
AhR-interacting congeners in the Aroclor 1254 mixture or due to the lack of Cyp expression in the brain. Very few studies on the presence of Cyp's in the brain have been reported and these generally report very low constitutive levels. Of the over 800 PubMed references on the Cyp1a2 gene, only 15 studies were in brain and only 3 of those investigated the AhR signaling pathway. The connection of the Col1a2 gene is via its regulation by TNF (Umannová et al., 2007), a well known AhR-regulated gene. Changes in expression of Col1a2, with its roles in cell migration, proliferation and adhesion, suggest an AhRmediated developmental effect. In regards to the ‘fifth’ criticism, it is indeed true that the same data set was analyzed differently in the two papers. However we fail to see how this is interpreted as a “flaw”. Genomic data sets are routinely analyzed and re-analyzed to pull out different types of information and these data are published separately. The difference here is that we chose to present the data in two companion manuscripts as the most efficient way to portray different aspects of this multi-dimensional study. The pathway analyses in the mode of action manuscript (Royland and Kodavanti, 2008) support and expand the developmental aspects outlined in the developmental manuscript (Royland et al., 2008). In the first manuscript (Royland et al., 2008) we evaluated global developmental effects by determining statistically significant changes in both constitutive gene expression and in gene expression following Aroclor 1254 treatment in two brain areas between PND7 and PND14. Genes in the individual condition lists are subsets of the total number of genes statistically different with age (see Figure 1, Venn diagram). In each case the top 100 differentially expressed genes were chosen as the most likely to be representative of the most impacted functional groups (complete gene lists are available from the journal on line as supplements). As it was our intent to examine functional correlates, we portrayed the relative number and direction of change for statistically different genes (PND7 to PND14) within functional groups for each condition and tissue (Figures 4 and 7). No statistics were necessary here as these lists are of genes identified that changed with age only in control or only in Aroclor-treated animals, but not in both. Likewise, the GO ontology analyses report p-values for each tissue/treatment and can be compared as to relative importance, but are not designed to report across treatment differences here. In any developmental study, developmental effects will be inseparable from treatment effects, thus it is not surprising that many genes were found in common between the two treatments. A similar functional analysis was carried out on these genes (Figures 3 and 6) to add greater depth to the study. Commentary on how the relative expression level was impacted for individual genes that were common between the conditions (Tables 2 and 4) was merely to highlight a possible general effect of Aroclor 1254 on gene expression. As stated in the figure legends and as is most often done in heatmap generation, the heatmaps display expression changes relative to the mean of all expression across all samples/conditions that are included in the heatmap (note, each figure contains heatmap representations of three individual clustering analyses). We believe all the information required to understand Figures 5 and 8 is included in the accompanying legends. We also believe that our analyses across developmental ages in this first manuscript (Royland et al., 2008), which is different from the more traditional comparison of control vs. treatment in the companion manuscript, is a novel approach that adds new insight into global Aroclor 1254 effects and its impact on neurodevelopment, providing possible new avenues for research. The concerns Drs. Carlson and Silkworth raised regarding control/ Aroclor comparisons were addressed in the companion manuscript on possible modes of action (Royland and Kodavanti, 2008). During the brain growth spurt occurring during much of the dosing regimen, individual brain areas are rapidly changing, each on its own timetable. This is apparent in the multivariate PCA analysis (Figure 9) showing that gene expression is variable by brain area, age and treatment (especially at PND7). These rapidly changing gene environments
259
resulted in an insufficient number of genes (see Venn diagrams, Figures 1 and 4) to do the kind of gene-by-gene pathway analysis across both time and dose as recommended by Drs. Carlson and Stillworth. However, we strongly believe that our grouping of transcripts into functional categories is appropriate. Individual gene expression is known to be affected by timing, homeostatic mechanisms, feedback loops, and a host of normal physiological processes. Reporting on the pathways that were enriched (i.e., with multiple transcripts involved) under our experimental conditions provides important information on possible modes of action. Despite the concerns raised by Drs. Carlson and Silkworth, we strongly believe that, although the current genomic data are preliminary, the results from this study are very useful in guiding complementary studies using proteomics and neurochemical measures to determine mode(s) of action for these chemicals. References Behar, T.N., Li, Y.X., Tran, H.T., Ma, W., Dunlap, V., Scott, C., Barker, J.L., 1996. GABA stimulates chemotaxis and chemokinesis of embryonic cortical neurons via calcium-dependent mechanisms. J. Neurosci. 16, 1808–1818. Behar, T.N., Scott, C.A., Greene, C.L., Wen, X., Smith, S.V., Maric, D., Liu, Q.Y., Colton, C.A., Barker, J.L., 1999. Glutamate acting at NMDA receptors stimulates embryonic cortical neuronal migration. J. Neurosci. 19, 4449–4461. Behar, T.N., Schaffner, A.E., Scott, C.A., Greene, C.L., Barker, J.L., 2000. GABA receptor antagonists modulate postmitotic cell migration in slice cultures of embryonic rat cortex. Cereb. Cortex 10, 899–909. Borlok, J., Dangers, M., Thum, T., 2002. Aroclor 1254 modulates gene expression of nuclear transcription factors: implications for albumin gene transcription and protein synthesis in rat hepatocyte cultures. Toxicol. Appl. Pharmacol. 181, 79–88. Bush, B., Kadlec, M.J., 1995. Dynamics of PCBs in the aquatic environment. Great Lakes Research Reviews 1, 24–30. Dong, H., Wade, M., Williams, A., Lee, A., Douglas, G.R., Yauk, C., 2005. Molecular insight into the effects of hypothyroidism on the developing cerebellum. Biochem. Biophys. Res. Commun. 330, 1182–1193. Eriksson, P., Fredriksson, A., 1996. Developmental neurotoxicity of four orthosubstituted polychlorinated biphenyls in the neonatal mouse. Environ. Toxicol. Pharmacol. 1, 155–165. Fiszman, M.L., Borodinsky, L.N., Neale, J.H., 1999. GABA induces proliferation of immature cerebellar granule cells grown in vitro. Brain Res. 115, 1–8. Hakak, Y., Walker, J.R., Li, C., Wong, W.H., Davis, K.L., Buxbaum, J.D., Haroutunian, V., Fienbert, A.A., 2001. Genome-wide expressin analysis reveals dysregulation of myelination-related genes in chronic schizophrenia. PNAS. 98, 4746–4751. Kaiser, S., Foltz, L.A., George, C.A., Kirkwood, S.C., Bemis, K.G., Lin, X., Gelbert, L.M., Nisenbaum, L.K., 2004. Phencyclidine-induced changes in rat cortical gene expression identified by microarray analysis: implications for schizophrenia. Neurobiol. Dis. 16, 220–235. Kodavanti, P.R.S., 2005. Neurotoxicity of persistent organic pollutants: possible mode(s) of action and further considerations. Dose Response 3, 273–305. Kodavanti, P.R.S., Kannan, N., Yamashita, N., Derr-Yellin, E.C., Ward, T.R., Burgin, D.E., Tilson, H.A., Birnbaum, L.S., 2001. Differential effects of two lots of Aroclor 1254: congener-specific analysis and neurochemical end points. Environ. Health Perspect. 109, 1153–1161. Lee, M.H., Kim, M., Lee, B.H., Kang, K.S., Kim, H.L., Yoon, B.I., Chung, H., Kong, G., Lee, M.O., 2008. Subchronic effects of valproic acid on gene expression profiles for lipid metabolism in mouse liver. Toxicol. Appl. Pharmacol. 226, 271–284. Mayes, B.A., McConnell, E.E., Neal, B.H., Brunner, M.J., Hamilton, S.B., Sullivan, T.M., Peters, A.C., Ryan, M.J., Toft, J.D., Singer, A.W., Brown, J.F., Jr., Menton, R.G., Moore, J.A., 1998. Comparative carcinogenicity in Sprague-Dawley rats of the polychlorinated biphenyl mixtures Aroclors 1016, 1242, 1254, and 1260. Toxicol. Sci. 41, 62–76. Mejia-Gervacio, S., Marty, A., 2006. Control of interneuron firing pattern by axonal autoreceptors in the juvenile rat cerebellum. J. Physiol. 571, 43–55. Royland, J.E., Kodavanti, P.R.S., 2008. Gene expression profiles following exposure to a developmental neurotoxicant, Aroclor 1254; pathway analysis for possible mode(s) of action. Toxicol. Appl. Pharmacol. 231, 179–196. Royland, J.E., Wu, J., Zawia, N.H., Kodavanti, P.R.S., 2008. Gene expression profiles in the cerebellum and hippocampus following exposure to a neurotoxicant, Aroclor 1254: developmental effects. Toxicol. Appl. Pharmacol. 231, 165–178. Safe, S., 1994. Polychlorinated biphenyls (PCBs): environmental impact, biochemical, and toxic responses, and implications for risk assessment. Crit. Rev. Toxicol. 24, 87–149. Schantz, S.L., Levin, E.D., Bowman, R.E., Heironimus, M.P., Laughlin, N.K., 1989. Effects of perinatal PCB exposure on discrimination-reversal learning in monkeys. Neurotoxicol. Teratol. 11, 243–250. Umannová, L., Zatloukalová, J., Machala, M., Krcmár, P., Májková, Z., Hennig, B., Kozubík, A., Vondrácek, J., 2007. Tumor necrosis factor-alpha modulates effects of aryl hydrocarbon receptor ligands on cell proliferation and expression of cytochrome P450 enzymes in rat liver “stem-like” cells. Toxciol. Sci. 99, 79–89. Wei, Y., Liu, Y., Wang, J., Tao, Y., Dai, J., 2008. Toxicogenomic analysis of the hepatic effects of perfluorooctanoic acid on rare minnows (Gobiocypris rarus). Toxicol. Appl. Pharmacol. 226, 285–297.
260
Response to Letter to the Editor
Wolff, M.S., 1985. Occupational exposure to polychlorinated biphenyls (PCBs). Environ. Health Perspect. 60, 133–138. World Health Organization, 1993. Environmental health criteria 140: polychlorinated biphenyls and terphenyls, In: Dobson, S., van Esch, G.J. (Eds.), nternational Programme on Chemical Safety, 2nd ed. World Health Organization, Geneva, Switzerland. Wu, Z., Irizarry, R.A., Gentleman, R., Murillo, F.M., Spencer, F., 2003. A Model Based Background Adjustment for Oligonucliotide Expression Arrays. Technical Report. John Hopkins University, Dept. Biostatistics Working Papers, Baltimore. Yang, J.-H., Derr-Yellin, E.C., Kodavanti, P.R.S., 2003. Alterations in brain protein kinase C isoforms following developmental exposure to a polychlorinated biphenyl mixture. Mol. Brain Res. 111, 123–135.
Joyce E. Royland Prasada Rao S. Kodavanti* U.S. Environmental Protection Agency, Neurotoxicology Division, B105-06, Research Triangle Park, NC 27711, USA E-mail address: [email protected] (P.R.S. Kodavanti). *Corresponding author. 27 January 2009