Patrícia Pereira

Hilda De Pablo

Carlos Vale

Mário Pacheco

Combined use of environmental data and biomarkers in fish () inhabiting a eutrophic and metal-contaminated coastal system -gills reflect environmental contamination

Keywords: Oxidative stress, fish, metals, coastal lagoon, temporal variability

An investigative biomonitoring study was carried out in a eutrophic coastal system with a moderate contamination by metals (Óbidos lagoon, Portugal), combining the evaluation of exposure concentrations with metals accumulation and oxidative stress responses in gills of the golden grey mullet (Liza aurata). Two contrasting seasons (winter and summer) were considered at three sites: Barrosa (BB) and Bom-Sucesso (BS) branches; Middle lagoon (ML).

Data on the water column pointed to a higher metals and nutrients availability at BB that was reflected in the higher metal concentrations in gills, particularly in winter. Similarly, oxidative stress responses demonstrated a pro-oxidant challenge at BB (winter and summer), which was corroborated by an integrated biomarker response index (IBR). Metal concentrations in gills were higher in summer than winter, reflecting the increased environmental concentrations in combination with elevated metabolic rates. Catalase (CAT), glutathione-S-transferase (GST), total glutathione (GSH t ) and lipid peroxidation (LPO) increases observed in winter at BB were related with metal accumulation, while summer enhancement of glutathione peroxidase (GPx), glutathione reductase (GR), GST and GSH t was associated with other stressors. Inter-site differences on the basis of IBR were more accentuated in winter. Gills can be considered as an important route of entry for contaminants and were demonstrated to reflect water contamination and are therefore useful in the context of environmental assessment.

INTRODUCTION

The suitability of fish gills for environmental biomonitoring was recommended by several recent works [START_REF] Fernandes | Oxidative stress response in gill and liver of Liza saliens, from the Esmoriz-Paramos Coastal Lagoon, Portugal[END_REF][START_REF] Dautremepuits | Gill and kidney antioxidant processes and innate immune system responses of yellow perch (Perca flavescens) exposed to different contaminants in the St. Lawrence River, Canada[END_REF][START_REF] Oliveira | Contamination assessment of a coastal lagoon (Ria de Aveiro, Portugal) using defence and damage biochemical indicators in gill of Liza aurata -An integrated biomarker approach[END_REF]. In fact, gills have a wide surface area that is continuously in contact with the external medium and are thus the main uptake route of contaminants from aqueous phase [START_REF] Playle | Modelling metal interactions at fish gills[END_REF]. On the other hand, fish gills are able to accumulate chemicals that were taken up by other exposure routes, due to their position between the venous and arterial circulation, thus receiving nearly all of the cardiac output [START_REF] Levine | CYP1A expression in liver and gill of rainbow trout following waterborne exposure: implications for biomarker determination[END_REF]. Moreover, gills play a key role in fish physiology, such as gaseous exchange, osmotic and ionic regulation, acid-base balance and nitrogenous waste [START_REF] Ahmad | DNA damage and lipid peroxidation vs. protection responses in the gill of Dicentrarchus labrax L. from a contaminated coastal lagoon (Ria de Aveiro, Portugal)[END_REF]. It is documented that fish gills respond earlier than kidney and liver under field contamination exposure [START_REF] Ahmad | Enzymatic and non-enzymatic antioxidants as an adaptation to phagocyte-induced damage in Anguilla anguilla L. following in situ harbor water exposure[END_REF][START_REF] Santos | Anguilla anguilla L. antioxidants responses to in situ bleached kraft pulp mill effluent outlet exposure[END_REF]. Its applicability is also enhanced by the fact that gills' epithelium suffers an increased cell division to counteract its exfoliation and erosion and by this feature gills usually reflect current exposures [START_REF] Mieiro | Mercury distribution in key tissues of fish (Liza aurata) inhabiting a contaminated estuary -implications for human and ecosystem health risk assessment[END_REF]. This seems to be particularly relevant when fish, with their non-sessile behaviour, are used as bioindicators.

Monitoring the state of coastal systems has traditionally been largely based on measurements of concentrations of contaminants in sediment and water column. However, this chemical monitoring alone is no longer considered meaningful or cost-effective. Additionally, it has become evident that measuring concentrations in the tissues does not provide the information needed to assess biological effects that contaminants may cause in organisms [START_REF] Lehtonen | Monitoring biological effects of pollution in the Baltic Sea: Neglected -but still wanted?[END_REF]. In fact, contaminants usually appear in the environment as very complex mixtures that can cause interactive effects on the biota impossible to evaluate only by means of chemical analyses. Thus, recent investigative monitoring studies have involved the combined measurement of environmental concentrations and body burdens, together with the evaluation of biomarkers, which represent early signs of biological effects [START_REF] Van Der Oost | Fish bioaccumulation and biomarkers in environmental risk assessment: a review[END_REF]. In this context, many pollutants are capable of inducing oxidative stress in aquatic animals by disturbing the antioxidants efficiency and enhancing the intracellular reactive oxygen species (ROS), which often prelude in DNA damage, lipid peroxidation (LPO) and enzyme inhibition [START_REF] Van Der Oost | Fish bioaccumulation and biomarkers in environmental risk assessment: a review[END_REF]. Antioxidant defence systems include antioxidant enzymes such as superoxide dismutase (SOD), catalase (CAT), glutathione peroxidase (GPx) and glutathione reductase (GR), as well as non-enzymatic antioxidants (e.g. glutathione -GSH) which have been extensively used as biomarkers of oxidative stress [START_REF] Ahmad | Anguilla anguilla L. oxidative stress biomarkers: an in situ study of freshwater wetland ecosystem (Pateira de Fermentelos, Portugal)[END_REF]2008;[START_REF] Guilherme | Antioxidant and biotransformation responses in Liza aurata under environmental mercury exposure -Relationship with mercury accumulation and implications for public health[END_REF]. LPO estimation has also been found to have a high predictive value as a biomarker of effect [START_REF] Guilherme | Antioxidant and biotransformation responses in Liza aurata under environmental mercury exposure -Relationship with mercury accumulation and implications for public health[END_REF]. To our knowledge, only a few studies have combined the measurement of oxidative stress in fish gills with concentrations of accumulated contaminants [START_REF] Hansen | Antioxidative stress proteins and their gene expression in brown trout (Salmo trutta) from three rivers with different heavy metal levels[END_REF][START_REF] Fernandes | Oxidative stress response in gill and liver of Liza saliens, from the Esmoriz-Paramos Coastal Lagoon, Portugal[END_REF].

Additionally, an investigative approach with the previous information and environmental contamination levels is more uncommon [START_REF] Hansen | Antioxidative stress proteins and their gene expression in brown trout (Salmo trutta) from three rivers with different heavy metal levels[END_REF].

Previous studies have employed the golden grey mullet (Liza aurata) as a model species on field studies [START_REF] Oliveira | Contamination assessment of a coastal lagoon (Ria de Aveiro, Portugal) using defence and damage biochemical indicators in gill of Liza aurata -An integrated biomarker approach[END_REF][START_REF] Guilherme | Antioxidant and biotransformation responses in Liza aurata under environmental mercury exposure -Relationship with mercury accumulation and implications for public health[END_REF] and under laboratory conditions [START_REF] Cionna | Effects of nonylphenol on juveniles and adults in the grey mullet, Liza aurata[END_REF][START_REF] Oliveira | Cytochrome P4501A, genotoxic and stress responses in golden grey mullet (Liza aurata) following short-term exposure to phenanthrene[END_REF]. In fact, L. aurata has favourable features as a sentinel of contaminated coastal systems. This is a pelagic species that is in frequent contact with sediments, filtering its superficial layer and suspended particulate matter (SPM) in the water column. Moreover, mullets play an important role in the estuarine trophic web [START_REF] Almeida | Feeding ecology of Liza ramada (Risso, 1810) (Pisces, Mugilidae) in a south-western estuary of Portugal[END_REF].

In order to access the environmental health of a coastal lagoon impacted by eutrophication and moderate metal contamination (Óbidos lagoon, Portugal), gills of Liza aurata were monitored for metal accumulation and oxidative stress in two contrasting seasons: winter and summer. Peroxidative damage and antioxidant (enzymatic and non-enzymatic) responses were measured in L. aurata gills and their relationship with metal accumulation as well as with environmental chemical data (metals and nutrients) was assessed. It was also intended to evaluate the suitability of the adopted combined strategy (external levels of exposure/bioaccumulation markers/oxidative stress responses) applied to fish gills on monitoring coastal ecosystems.

MATERIALS AND METHODS

Study area

The Óbidos lagoon is a shallow coastal ecosystem, located on the west coast of Portugal with a wet area of 7 km 2 , permanently connected to the sea through a narrow inlet (Figure 1). It comprises areas of different morphological and sedimentary characteristics: sand banks and narrow channels in the lower/middle lagoon; muddy bottom sediments in the two inner branches (Barrosa and Bom-Sucesso). The Barrosa branch is shallower (mean depth 0.5-1 m) compared to the Bom-Sucesso (mean depth 1-2 m) and water circulation is mostly driven by tides and by a small tributary (Cal River) that drains agricultural fields. Urban effluents from a nearby town (Caldas da Rainha, 50 000 inhabitants) also have been discharged to the Barrosa branch by the Cal River. Consequently, this area presents the highest nutrient availability of the lagoon, being classified as eutrophic (Pereira et al., 2009a). High nutrient concentrations were in line with abundant macroalgae (Ulva sp. and Enteromorpha sp.) and a broad daily variation of dissolved oxygen concentration during the summer months [START_REF] Pereira | Metal and nutrient dynamics in a eutrophic coastal lagoon (Óbidos, Portugal): the importance of observations at different time scales[END_REF]. Previous water quality surveys of the lagoon showed maximum concentrations of Mn, Ni, Cu and Cd in Barrosa branch. Major metal sources were identified to be the Cal River in periods of higher inflow and remobilization from sediments in summer months [START_REF] Pereira | Metal and nutrient dynamics in a eutrophic coastal lagoon (Óbidos, Portugal): the importance of observations at different time scales[END_REF]. The Bom-Sucesso branch is also a confined area but receives a smaller freshwater flow (Vala do Ameal) with better water quality than the Cal River, according to the Portuguese categorization of freshwater systems. In this area, metal remobilization from sediments is less likely due to the greater depth.

The middle lagoon was previously characterised by coarser sediments with a low affinity for metals [START_REF] Carvalho | The use of the marine biotic index AMBI in the assessment of the ecological status of the Óbidos lagoon (Portugal)[END_REF], as well as by a better water quality (Pereira et al., 2009a).

Moreover, previous ecotoxicological studies with Ulva sp. and Carcinus maenas pointed out Barrosa branch as a critical area of the lagoon with potential impact on autochthonous populations, whereas the lower/middle lagoon was considered the reference area (Pereira et al., 2009c;2009d).

Sampling

Two surveys were carried out at Óbidos lagoon in winter (December 2006) and summer (July 2007) during low-tide and juveniles of the golden grey mullet (L. aurata) were collected (n=10) using a traditional beach-seine net named "chincha". Three sampling areas were selected taking into account previous ecotoxicological studies (Pereira et al., 2009c(Pereira et al., , 2009d) ) (Figure 1): Middle lagoon (ML); Bom-Sucesso (BS) and Barrosa (BB) located in the inner branches of the lagoon. Fish biometrical parameters, such as weight and length ranged from 20-40 g and 13-17 cm. Juveniles were selected in order to minimize any gender interference and reduce the background arising from previous exposures to contaminants [START_REF] Ahmad | DNA damage and lipid peroxidation vs. protection responses in the gill of Dicentrarchus labrax L. from a contaminated coastal lagoon (Ria de Aveiro, Portugal)[END_REF]. Immediately after catching, fish were sacrificed and gills were removed and frozen in liquid nitrogen in two sets of samples: one for oxidative stress quantifications and other corresponding to the rest of the organ for metal determinations. In the laboratory, samples were preserved until further processing at -80 and -20 °C for oxidative stress and metal determinations, respectively.

Surface water (0.2 m depth) was sampled in the three areas for determinations of nutrients, chlorophyll a and metals.

Determination of physicochemical parameters and metals in water

Water temperature, salinity and dissolved oxygen were measured in situ using a YSI, 650 meter. Water samples for determinations of nitrate+nitrite (NO 3 -+NO 2 -), ammonium (NH 4 + ), phosphate (PO 4 3-) and silicate Si(OH) 4 ) were filtered through MSI Acetate Plus filters and analyses carried out using an autoanalyser TRAACS 2000 (Bran+Luebbe). For chlorophyll a determinations, water was filtered through a Whatman GF/F (0.7 µm) filter that was immediately frozen at -20 ºC and later extracted in 90% acetone, for analysis in a Perkin Elmer Fluorometer using the protocol modified by [START_REF] Lorenzen | A method for the continuous measurement of in vivo chlorophyll concentration[END_REF].

Manganese, Cu, Ni, Cr, Pb and Cd in the collected waters (triplicate samples) were measured using diffusive gradients of thin films (DGT). All DGT holders, Chelex-100 resins and diffusive gels (type APA, 0.8-mm thickness, open pore> 5 nm) [START_REF] Zhang | Diffusional characteristics of hydrogels used in DGT and DET techniques[END_REF] were purchased from DGT Research (Lancaster, UK). The DGT devices were deployed in 2-L polypropylene bottles with unfiltered sampled water and stirred at constant temperature for 48 hours. After devices retrieval, resins were eluted by immersion in 5 mL of 1 M HNO 3 (prepared from suprapur nitric acid) at a minimum of 24 h. Eluates were analysed directly by a quadropole inductively coupled plasma mass spectrometer (ICP-MS) (Thermo Elemental, X-Series). All eluates were analysed with reagents blanks and an international standard of river water (SLRS-4) used to control eventual contaminations during the analytical procedure and the procedure accuracy, respectively. The eluate concentration was converted into the mass of metal accumulated on the resin using an elution factor with a yield value of 0.8 [START_REF] Zhang | Performance-characteristics of diffusion gradients in thinfilm for the in situ measurement of trace-metal in aqueous-solution[END_REF] and a resin gel volume of 0.15 mL. DGT-metal concentrations were calculated according to the following equation: [Metal]=(M g)/(tAD); where M is the mass of metal accumulated on the resin during the emersion time (t); ∆g the diffusive gel thickness (0.08 cm); A the exposure area (3.14 cm 2 ); and D the diffusion coefficient of the metal in the gels as provided by DGT Research.

Water was collected in triplicates for separation of suspended particulate matter (SPM), which was retained on 0.45 µm membranes that were mineralized completely with HF (40%)

and Aqua Regia (HCl-36%:HNO 3 -60%; 3:1) in closed Teflon bombs (100ºC for 1 h), evaporated to near dryness (DigiPrep HotBlock -SCP Science), redissolved with 1 mL of doubled-distilled HNO 3 and 5 mL of ultra-pure water, heated for 20 min at 75 ºC, and diluted to 50 mL with ultrapure water. Reagents blanks and international certified standards of sediments from the National Research Council of Canada (1646a; BCSS-1; MESS-3) were prepared in a similar way of samples to control the accuracy of the procedure. The concentrations of Mn and Zn were determined by flame atomic absorption spectrometry (Perkin Elmer Analyst 100) whereas Cu, Ni, Cr, Pb and Cd were determined by ICP-MS. Levels of the analysed elements obtained in the reference materials were consistently within the ranges of certified values.

Metal analyses in gills

Approximately 50 mg of freeze dried tissue was digested with a mixture of HNO 3 (doubled distilled from 65%) and H 2 O 2 (suprapure, 30%) at 60 °C for 12 h, at 100 °C for 1 h and at 80 ºC for 1 h according to the method described in Pereira et al. (2009e). Procedural blanks were prepared using the same analytical procedure and reagents. Concentrations of Mn, Zn, Cu, Ni, Cr, Pb and Cd were determined by ICP-MS. International certified standards (TORT-2, DOLT-2 and DOLT-3) were use to control the accuracy of the analytical procedure.

Biochemical analyses in gills

Tissue handling and preparation of post-mitochondrial supernatant (PMS) fraction -Tissue samples were homogenized, using a Potter-Elvehjem homogenizer, in chilled phosphate buffer (0.1 M, pH 7.4) (1 g of tissue/10 ml buffer). This homogenate was then divided in two aliquots for LPO and PMS preparation. The PMS preparation was obtained by centrifugation in a refrigerated centrifuge (Optima TL, Beckman) at 13,400 g for 20 min at 4 ºC.

Aliquots of PMS were divided in microtubes and stored at -80 ºC until analyses.

Estimation of lipid peroxidation (LPO) -LPO was determined in the previously prepared homogenate as adapted by [START_REF] Filho | Seasonal changes in antioxidant defences of the digestive gland of the brown mussel (Perna perna)[END_REF]. Briefly, to 150 µl of homogenate, 10 µl of 1-1 butylated hydroxytoluene (4% in methanol) was added and mixed well to prevent oxidation. The absorbance of each sample was measured at 535 nm. The rate of LPO was expressed in nmol of thiobarbituric acid reactive substances (TBARS) formed per milligram of fresh tissue using a molar extinction coefficient of 1.56x10 5 M -1 cm -1 .

Catalase (CAT) activity measurement -CAT activity was assayed in PMS by

Claiborne ( 1985) method (at 25 ºC). Briefly, the assay mixture consisted of 1.95 ml phosphate buffer (0.05 M, pH 7.0), 1 ml hydrogen peroxide (0.019 M) and 50 µl of sample in final volume of 3 ml. Change in absorbance was recorded at 240 nm and CAT activity was calculated in terms of µmol H 2 O 2 consumed min -1 mg -1 protein using a molar extinction coefficient of 43.5 M -1 cm -1 .

Glutathione peroxidase (GPx) activity measurement -GPx activity was determined in PMS according to the method described by [START_REF] Mohandas | Differential distribution of glutathione and glutathione related enzymes in rabbit kidney. Possible implications in analgesic neuropathy[END_REF] with some modifications (at 25 ºC). The assay mixture consisted of 0.72 ml phosphate buffer (0.05 M, pH 7.0), 0.05 ml EDTA (1 mM), 0.05 ml sodium azide (1 mM), 0.025 ml GR (1 IU/ml), 0.05 ml GSH (4 mM), 0.05 ml NADPH (0.8 mM), 0.005 ml H 2 O 2 (1.0 mM) and 0.05 ml of PMS in a total volume of 2 ml.

Oxidation of NADPH was recorded at 340 nm, and GPx activity was calculated in terms of nmol NADPH oxidized min -1 mg -1 protein using a molar extinction coefficient of 6.22x10 3 M -1 cm -1 .

Glutathione reductase (GR) activity measurement -GR activity was assayed by the method of [START_REF] Cribb | Use of a microplate reader in an assay of glutathione reductase using 5,5'-dithiobis(2-nitrobenzoic acid)[END_REF] with some modifications (at 25 ºC). Briefly, the assay mixture contained 0.025 ml of PMS fraction and 0.925 ml of NADPH (0.2 mM), GSSG (1 mM) and DTPA (0.5 mM). The enzyme activity was quantified by measuring the disappearance of NADPH at 340 nm during 3 min. The enzyme activity was calculated as nmol NADPH oxidized min -1 mg -1 protein using a molar extinction coefficient of 6.22×10 3 M -1 cm -1 .

Glutathione-S-transferase (GST) activity measurement -GST activity was determined using CDNB (1-chloro-2,4-dinitrobenzene) as substrate according to the method of [START_REF] Habig | Glutathione S-transferases. The first enzymatic step in mercapturic acid formation[END_REF] (at 25 ºC). The assay was carried out in a quartz cuvette with a 2 ml mixture of phosphate buffer (0.2 M, pH 7.4), CDNB (0.2 mM) and 0.2 M GSH. The reaction was initiated by addition of 10 µl of PMS and the increase in absorbance was recorded at 340 nm during 3 min. The enzyme activity was calculated as nmol CDNB conjugate formed min -1 mg -1 protein using a molar extinction coefficient of 9.6 mM -1 cm -1 .

Total glutathione (GSH t ) content measurement -Protein content in the tissue homogenate was precipitated with TCA (5%) for 1 h and then centrifuged at 13,400 g for 20 min at 4 ºC. GSH t content was determined in the resulting supernatant (deproteinated PMS) (at 25 ºC) adopting the enzymatic recycling method using glutathione reductase (GR) excess, whereby the sulfhydryl group of GSH reacts with DTNB (Ellman's reagent) producing a yellow coloured 5-thio-2-nitrobenzoic acid (TNB). The rate of TNB production is directly proportional to this recycling reaction, which is in turn directly proportional to the concentration of GSH in the sample. Formation of TNB was measured by spectrophotometry at 412 nm. It should be noted that GSSG is converted to GSH by GR in this system, which consequently measures total GSH.

The results were expressed as nmol TNB formed min -1 mg -1 protein using a molar extinction coefficient of 14.1×103 M -1 cm -1 .

Protein measurement -Total protein contents were determined according to the Biuret method [START_REF] Gornall | Determination of serum proteins by means of the biuret reaction[END_REF], using bovine serum albumin (E. Merck-Darmstadt, Germany) as a standard.

Statistical analyses

A method for combining all the measured oxidative stress endpoints into a general index -Integrated Biomarker Response (IBR) [START_REF] Beliaeff | Integrated biomarker response: a useful tool for ecological risk assessment[END_REF] -was applied to data of both surveys. The basis of data processing is described here briefly. First, the mean value (X i ) for each biomarker at each sampling area and season was calculated. In addition, the general mean (m i ) and standard deviations (s i ) of each biomarker were estimated for all sampling areas and seasons. The value of X i was then standardized to obtain Y i , where Y i = (X im i )/s i . Then Z i was computed via the equation Z i = -Y i or Z i = Y i in the case of a biological effect corresponding, respectively, to inhibition or activation. The minimum value (Min i ) of Z i for each biomarker was calculated for all sampling areas and seasons, and then the score S i was computed as S i = Z i +Min i , where Min i  is the absolute value. Finally, IBR for each sampling site and season was calculated via de following formula: IBR= (S 1 xS 2 )/2 + (S 2 xS 3 )/2 +…+ (S n- 1 xS n )/2 + (S n xS 1 )/2; in which the obtained score for each biomarker (S i ) is multiplied with the score of the next biomarker (S i+1 ), arranged as a set, dividing each calculation by 2 and summing-up of all values.

Statistical software (Statistica 7.0) was used for statistical analyses. The assumptions of normality and homogeneity of data were verified. ANOVA analysis was used to compare sampling sites and Tukey test was applied for Post-hoc comparison [START_REF] Zar | Biostatistical Analysis[END_REF]. Differences between means were considered significant when p<0.05. The correlation between metals in the water column (DGT-measured and SPM fractions) and accumulated in gills, as well as between accumulated metals and oxidative stress responses in gills was tested (considering the mean values of the two surveys at each site) by Pearson analysis (p<0.05).

RESULTS

Water characteristics

Water temperature was slightly higher at BB and BS in summer and winter, while salinity was lower at BB than at the other sites, particularly in winter (Table 1). In general, dissolved oxygen was around 100% but undersaturation was recorded at BB in summer and BS in winter. In general, the maximum values of nitrogenous compounds, Si(OH) 4 , PO 4 3-and chlorophyll a were recorded at BB in both surveys. NH 4 + , NO 3 -+NO 2 -and Si(OH) 4 exhibited a typical temporal fluctuation in the lagoon with lower values in summer. These site and seasonal differences were consistent with previous measurements (Pereira et al., 2009a).

Table 2 presents the levels of metals in water (DGT-measured fraction) at the three sites in winter and summer. BB presented significantly higher levels of Cu in winter than ML. In summer, Cu was significantly higher at BB (than ML and BS), while the opposite pattern was observed for Mn. Moreover, in summer Cr was significantly higher at BB and BS than at ML. Concerning temporal differences, significantly higher levels of Mn and Cd were found in winter, while Cu and Cr enhanced in summer.

Copper and Cd in SPM collected in winter were significantly higher at BB (than ML, ML and BS, respectively) (Table 3). In summer, Cu, Ni and Pb in SPM were higher at BB, while SPM from ML contained higher levels of Mn, Zn and Cd than other sites (though without statistical significance). Moreover, Ni increased significantly in summer in relation to winter (BB)

whereas Pb (ML) and Cd (BB) were significantly higher in winter than in summer.

Metals in L. aurata gills

In winter, gills of fish from BB showed significantly higher accumulation of Zn, Cu, Cr and Cd than ML (Figure 2). Moreover, concentrations of Cr were significantly higher in gills of fish from BB than BS in winter, whereas at the BS site Zn and Pb were elevated in comparison with ML. In summer, gills at BB showed lower concentrations of Cu (than BS) and Cr (than BS and ML).

Differences between the two surveys were found for all metals, except Mn and Pb. Gills in winter presented higher contents of Cu and Ni (BB) than in summer. In contrast, in summer gills presented higher concentrations of Zn (all sites), Cu (BS), as well as Cr and Cd (ML, BS).

Oxidative stress in L. aurata gills

A spatial analysis of gills antioxidant defences revealed significant increases of CAT (winter), GPx and GR (summer) activities, as well as GSH t content (winter and summer) at BB relatively to ML and BS (Figure 3). Moreover, GST activity in winter and summer was enhanced significantly at BB in comparison to ML. In terms of peroxidative damage, significantly higher levels were also found at BB (versus ML and BS) in winter. Concerning temporal differences, significantly higher levels were found in summer for: CAT, GPx and GR at all sites; GST at ML and BS; and GSH t at BB. In addition, LPO was higher in winter than summer at BB.

Correlations between parameters and IBR calculation

The significant correlations obtained between metal levels in gills versus oxidative stress endpoints are shown in Figure 4 (non-significant correlations were not presented). Positive correlations were found between activities of GPx, GR and GST and accumulated levels of Zn, Cd and Cu, respectively. Significant correlations were not found between metals in the water column (DGT-measured fraction and SPM) and accumulated in L. aurata gills.

The antioxidants and LPO values measured in gills were processed in order to calculate the IBR index (Figure 5). IBR values showed different variation ranges in winter and summer, being higher in summer. Even so, in both surveys the computed IBR values were always higher at BB than other sites. In particular, the estimated IBR value for BB was 100 times higher than ML in winter, whereas in summer the difference between the two sites was approximately 3 times.

DISCUSSION

Water quality and metal availability

In winter and summer BB exhibited a higher metal availability than BS and particularly ML. This was particular evident for Cu in water (DGT-measured fraction) and SPM, Cr in water, as well as Mn, Ni and Cd in SPM. The tributary that discharges into the BB branch probably contributed to metal inputs, mainly in periods of higher freshwater flows [START_REF] Pereira | Metal and nutrient dynamics in a eutrophic coastal lagoon (Óbidos, Portugal): the importance of observations at different time scales[END_REF].

In fact, levels of water Cd and Mn (all sites) as well as particulate Pb (all sites) and Cd (BB)

were higher in winter than in summer. Conversely, values of Cu (all sites) and Cr (BS, BB) in water were higher in summer than winter, when freshwater flows into the lagoon are almost insignificant. Nickel associated with particles also increased in summer at BB. These summer enhancements are probably related with the pulse release of metals from sediments during periods of low oxygenation as indicated by the simultaneous undersaturation of oxygen recorded mainly at BB. This internal source was hypothesized to occur in BB branch [START_REF] Pereira | Metal and nutrient dynamics in a eutrophic coastal lagoon (Óbidos, Portugal): the importance of observations at different time scales[END_REF]. The importance of these sediment-related processes for the water quality is probably lower at BS and ML than at BB due to the higher trophic state of these two areas [START_REF] Pereira | Metal and nutrient dynamics in a eutrophic coastal lagoon (Óbidos, Portugal): the importance of observations at different time scales[END_REF]. In fact, BB was already proposed as a major impacted area in the lagoon being classified as eutrophic (Pereira et al., 2009a). Accordingly, BB exhibited the highest levels of nutrients and chlorophyll a in both seasons. In wet seasons, the tributary that discharges at BB branch introduce large amounts of nutrients into that area of the lagoon (Pereira et al., 2009a). Since PO 4 3-regeneration could occur in organically-rich sediments with upwards diffusion during periods of low oxygenation (mainly in summer), this should be considered as an additional input to the water column at BB as previously recorded [START_REF] Pereira | Metal and nutrient dynamics in a eutrophic coastal lagoon (Óbidos, Portugal): the importance of observations at different time scales[END_REF].

Metal concentrations in L. aurata gills and associations with environmental factors

Gills are the first organ to be in contact with water and resuspended sediment particles, so they can be relevant sites of interaction with metal ions [START_REF] Fernandes | Bioaccumulation of heavy metals in Liza saliens from the Esmoriz-Paramos coastal lagoon, Portugal[END_REF]. In winter, gills from BB presented higher concentrations of Cu, Cr and Cd in comparison with ML, which indicates a higher uptake from water column and is in line with the spatial patterns recorded for those metals in SPM and water (DGT-measured fraction). This gills feature was previously observed in other fish species [START_REF] Reynders | Accumulation and effects of metals in caged carp and resident roach along a metal pollution gradient[END_REF]. The Pb enhancement at BS in winter (in comparison with ML) is also probably a response to a higher availability as suggested by levels in SPM.

Gills from BB and BS showed higher contents of Zn in winter but this was not accompanied by increases in concentrations in the water or SPM, which suggests other factors explaining its accumulation rather than simply Zn availability. This discrepancy was already observed by other authors for essential elements in fish under field exposure [START_REF] Demirak | Heavy metals in water, sediment and tissues of Leuciscus cephalus from a stream in southwestern Turkey[END_REF][START_REF] Reynders | Accumulation and effects of metals in caged carp and resident roach along a metal pollution gradient[END_REF]. An inconsistency between environmental data and accumulation was also found in summer for Cu and Cr. Despite the higher availability of those elements at BB their uptake and/or accumulation was lower than at other sites.

Manganese accumulation in gills did not reflect the spatial differences recorded in water or particles in summer. These results suggested that L. aurata was able to regulate tissue Mn concentrations as previously reported in the Óbidos lagoon for crabs (Pereira et al., submitted).

Indeed, the concentrations of essential metals (as Mn) in organisms tend to be highly regulated compared to nonessential. Fish, in particular, can use distinct strategies of metal homeostasis to achieve a steady-state balance: reducing metal accumulation and toxicity including uptake inhibition, increased elimination and detoxification, and storage [START_REF] Fernandes | Bioaccumulation of heavy metals in Liza saliens from the Esmoriz-Paramos coastal lagoon, Portugal[END_REF].

Pearson analysis did not provide any association between accumulated metals in gills and the corresponding levels in water or particles. Since data of the two seasons were considered jointly in the Pearson analysis and the associations above discussed did not occur in summer, non-significant correlations were obtained. The absence of significant correlations between water data and concentrations in the gills may have been due to the discrete water sampling compared to the time-integrated measurements in organisms.

In general, gills exhibited higher metal concentrations in summer. For instance, Cr was higher in summer (ML and BS) in line with the increases recorded in water (DGT-measured fraction). Identically, a higher availability could also be on the basis of Cu enhancement in gills during the summer at BS. However, gills concentrations of Zn and Cd in the gills increased in summer in contrast to relatively constant concentrations in the environment or enhancements in Cd in the winter. Considering that Zn is an essential element [START_REF] Fernandes | Bioaccumulation of heavy metals in Liza saliens from the Esmoriz-Paramos coastal lagoon, Portugal[END_REF], the enhancement of its uptake and subsequent accumulation in gills could represent an adaptation strategy in order to compensate for a higher requirement during the summer period.

Nonetheless, it is well known that the rates of metal uptake and accumulation increase with increasing temperature in ectothermal organisms [START_REF] Sokolova | Interactive effects of metal pollution and temperature on metabolism in aquatic ectotherms: implications of global climate change[END_REF]). An increase in metabolic rates at elevated temperatures may contribute to metal accumulation in ectotherms due to a higher energy demand, which results in elevated ventilation [START_REF] Pörtner | Climate variations and the physiological basis of temperature dependent biogeography: systemic to molecular hierarchy of thermal tolerance in animals[END_REF]. Under conditions of moderate contamination, as in the Óbidos lagoon, the seasonal difference of metabolism probably dominated and masked any differences in metal environmental availability, thus explaining the limited inter-site differences in gill metal concentrations observed in the summer.

Oxidative stress in L. aurata gills and causative factors

Spatial variations

Enzymatic antioxidant responses showed that fish are experiencing from pro-oxidant challenges at BB in winter and summer. This site presents compounds that stimulate the production of H 2 O 2 as depicted by the increase of CAT (winter) and GPx activities (summer).

Hydrogen peroxide is the main cellular precursor of the hydroxyl radical that is the most toxic ROS. Those enzymes are known to protect the cell by reducing H 2 O 2 to H 2 O. The induction of both enzymes has been described previously in gills of L. aurata exposed to contaminants under field conditions [START_REF] Oliveira | Contamination assessment of a coastal lagoon (Ria de Aveiro, Portugal) using defence and damage biochemical indicators in gill of Liza aurata -An integrated biomarker approach[END_REF]. In particular, metals were previously associated with CAT induction in fish gills [START_REF] Fernandes | Oxidative stress response in gill and liver of Liza saliens, from the Esmoriz-Paramos Coastal Lagoon, Portugal[END_REF]. In the Óbidos lagoon, CAT enhancement at BB in winter may be associated with higher accumulation of metals (particularly Cu, Cr and Cd). On the other hand, it is unlikely that GPx induction at BB in summer could be due to accumulated metals since substantial increases were not recorded in gills. In this direction, the significant correlation found between GPx and Zn accumulated in gills (provided by data of the two seasons jointly) has a limited value in establishing causal relationships strictly in summer. It is more plausible that other stressors rather than metals could be the basis of GPx induction in summer. The broad fluctuations of dissolved oxygen between day (180%) and night (30%) that were recorded in summer (Pereira et al., 2009a) could lead to an activation of oxidative stress defences. In fact, it is well known that organisms exposed to fluctuating levels of oxygen are subjected to marked oxyradical burst, which can be counteracted by elevated levels of antioxidants to prevent oxidative damage [START_REF] Storey | Oxidative stress: animal adaptation in nature[END_REF].

GR plays a major role in GPx and GST reactions as an adjunct in the control of peroxides and free radicals by maintaining the GSH redox status [START_REF] Bompart | Rapid automated analysis of glutathione reductase, peroxidase, and S-transferase activity: application to cisplatin-induced toxicity[END_REF]. Thus, exposure to pro-oxidant stressors leads to elevated activities of GR in aquatic organisms [START_REF] Regoli | Application of biomarkers for assessing the biological impact of dredged materials in the Mediterranean: the relationship between antioxidant responses and susceptibility to oxidative stress in the red mullet (Mullus barbatus)[END_REF] as in L. aurata gills [START_REF] Oliveira | Contamination assessment of a coastal lagoon (Ria de Aveiro, Portugal) using defence and damage biochemical indicators in gill of Liza aurata -An integrated biomarker approach[END_REF]. This is also occurring at BB site of the Óbidos lagoon in summer but was not associated with concomitant elevation of metals in the gills. The high GR activity at BB reflects increased glutathione recycling (replenishing GSH in order to avoid its depletion), suggesting that the ratio GSSG/GSH has been increased as consequence of pro-oxidants. As stated for GPx induction at BB in summer, other stressors rather than metals are probably on the basis of GR response. Nevertheless, GR was correlated with the Cd accumulated in the gills. As previously mentioned for GPx, this correlation has limited value in explaining summer responses.

In response to environmental changes, GSTs (a multigenic superfamily of multifunctional enzymes) may play a dual protective role associated with their activity on conjugation of electrophilic compounds (or phase I metabolites) with GSH ( [START_REF] Van Der Oost | Fish bioaccumulation and biomarkers in environmental risk assessment: a review[END_REF] and to a direct antioxidant action carried out by GST -class catalyzing the reduction of organic hydroperoxides by GSH [START_REF] Wang | Endogenous Glutathione Conjugates: Occurrence and Biological Functions[END_REF]. Given that the adopted methodology (using CDNB, which is conjugated by all GST isoforms with the exception of the q-class) determines total GST activity, the observed GST activity increase at BB in summer, compared with other sites, is not easily interpreted. When examining significant correlations (Pearson analysis), the GST induction was correlated with levels of Cu in the gills. In fact, several studies have shown that exposure to metals can lead to an increase of GST in liver [START_REF] Amado | Biomarkers in croakers Micropogonias furnieri (Teleostei: Sciaenidae) from polluted and non-polluted areas from the Patos Lagoon estuary (Southern Brazil): Evidences of genotoxic and immunological effects[END_REF][START_REF] Guilherme | Antioxidant and biotransformation responses in Liza aurata under environmental mercury exposure -Relationship with mercury accumulation and implications for public health[END_REF][START_REF] Oliveira | European eel (Anguilla anguilla L.) metallothionein, endocrine, metabolic and genotoxic responses to copper exposure[END_REF] and to a lesser extent in gills [START_REF] Dautremepuits | Gill and kidney antioxidant processes and innate immune system responses of yellow perch (Perca flavescens) exposed to different contaminants in the St. Lawrence River, Canada[END_REF].

GSH act against ROS protecting the cells from the pernicious effects of metals and other electrophilic compounds [START_REF] Wang | Endogenous Glutathione Conjugates: Occurrence and Biological Functions[END_REF]. In this study, GSH t was higher at BB in winter and summer, which indicates the presence of contaminants challenging cell defences depending on this thiol. Indeed, the variation of the cellular glutathione status in fish is considered an indicator of the degree and duration of exposure to oxidant pollutants [START_REF] Dautremepuits | Gill and kidney antioxidant processes and innate immune system responses of yellow perch (Perca flavescens) exposed to different contaminants in the St. Lawrence River, Canada[END_REF]. The protective and adaptive role of GSH was previously observed in L. aurata gills from a contaminated ecosystem [START_REF] Oliveira | Contamination assessment of a coastal lagoon (Ria de Aveiro, Portugal) using defence and damage biochemical indicators in gill of Liza aurata -An integrated biomarker approach[END_REF]. Although the Pearson analysis did not provide any significant relationships between accumulated metals and GSH t , the interference of metals (mainly Cu, Cr and Cd) in winter should be considered. In fact, metal metabolism involves the formation of GSH-metal complexes [START_REF] Atli | Responses of metallothionein and reduced glutathione in a freshwater fish Oreochromis niloticus following metal exposures[END_REF]. However, in summer other compounds are probably challenging cell defences depending on this thiol.

Elevated GSH t levels represent a great ability to destroy ROS and are a result of a transfer from other organs [START_REF] Oliveira | European eel (Anguilla anguilla L.) metallothionein, endocrine, metabolic and genotoxic responses to copper exposure[END_REF]. Indeed, most of the GSH is synthesized in liver and then exported to other tissues [START_REF] Deneke | Regulation of cellular glutathione[END_REF]. In addition, the concomitant increase of GR and GSH t (summer) indicates that the turnover of reduced glutathione pools was not enough to maintain the redox cycling and cell homeostasis, forcing a translocation from other organs to the gills. Moreover, this evidence of high GSH requirements is in agreement with the observed concomitant increases in GPx and GST activities.

Free radical reactions in biological membranes could form lipid hydroperoxides that decompose double bonds of unsaturated fatty acids, damaging lipid membranes ( [START_REF] Van Der Oost | Fish bioaccumulation and biomarkers in environmental risk assessment: a review[END_REF]. This complex process is known by lipid peroxidation. Our study revealed higher LPO at BB in winter, which is an indication of oxidative stress. This suggests that the observed mobilization of antioxidant defences was not able to prevent the gill LPO, probably induced by metal exposure in winter. Nevertheless, no associations were established between LPO and accumulated metals.

The higher levels of NH 4 + , NO 3 -+ NH 2 -and PO 4 3-registered at BB in winter could be associated with the observed oxidative stress responses. In fact, it has been suggested that ammonia and nitrite can be important pro-oxidants in fish, mussels and aquatic plants, leading to nitric radicals' production [START_REF] Lima | Biochemical responses of the marine mussel Mytilus galloprovincialis to petrochemical environmental contamination along the North-western coast of Portugal[END_REF][START_REF] Nimptsch | Ammonia triggers the promotion of oxidative stress in aquatic macrophyte Myriophyllum mattogrossense[END_REF]. At the Óbidos lagoon, the importance of those compounds on oxidative stress induction in crabs was previously invoked (Pereira et al., 2009d).

By combining the different biomarker signals, the IBR index provides a simple mean for a general description of populations' health status [START_REF] Broeg | Indices for the assessment of environmental pollution of the Baltic Sea coasts: Integrated assessment of a multi-biomarker approach[END_REF]. In the present work, IBR clearly demonstrated and confirmed that BB is the most impacted site within the lagoon.

Seasonal variations

In general, oxidative stress responses were enhanced in summer, namely activities of CAT, GPx, GR and GST and GSH t content. The sole exception to this temporal pattern occurred for LPO at BB. It should be also emphasised that despite the higher values in summer, GST and GSH t were the only endpoints exhibiting and identical spatial patterns in both seasons (characterised by increases at BB).

Metal concentrations in tissues are likely to change with the season, reflecting variability in the environmental inputs, but also changes in the metabolism [START_REF] Regoli | Trace Metals and Antioxidant Enzymes in Gills and Digestive Gland of the Mediterranean Mussel Mytilus galloprovincialis[END_REF]. Herein, seasonal variations of metals in gills were mainly associated with the effect of temperature on metals uptake, which enhanced in summer. Maximum enzymatic activities and GSH t content were recorded in summer when higher concentrations of Zn (at all sites) and Cd (at ML and BS)

were found. These results point to a major influence of metals for the recorded temporal variations.

Nonetheless, since water temperature rose 11 to 13 ºC in summer in relation to winter, its influence on oxidative stress responses is a possibility. Previous work carried out at the Óbidos lagoon, with crabs, also demonstrated activation of oxidative stress responses in summer (Pereira et al., submitted). Low temperature reduces the metabolic rates in ectothermal organisms, and hence lower enzymatic activities are in general observed in colder seasons.

This was reported for CAT and GST in fish [START_REF] Gorbi | Seasonal Variability of Metallothioneins, Cytochrome P450, Bile Metabolites and Oxyradical Metabolism in the European Eel Anguilla anguilla L. (Anguillidae) and Striped Mullet Mugil cephalus L. (Mugilidae)[END_REF][START_REF] Amado | Biomarkers in croakers Micropogonias furnieri (Teleostei: Sciaenidae) from polluted and non-polluted areas from the Patos Lagoon estuary (Southern Brazil): Evidences of genotoxic and immunological effects[END_REF]. Moreover, higher temperature lead to an increase in oxygen consumption and consequently to an enhancement of ROS generation [START_REF] Amado | Biomarkers in croakers Micropogonias furnieri (Teleostei: Sciaenidae) from polluted and non-polluted areas from the Patos Lagoon estuary (Southern Brazil): Evidences of genotoxic and immunological effects[END_REF]. Most aquatic organisms are able to adjust antioxidant defences as a mechanism to temperature adaptation [START_REF] Filho | Seasonal changes in antioxidant defences of the digestive gland of the brown mussel (Perna perna)[END_REF]. The current GSH t increase in summer is in agreement with previous studies with ectothermal marine animals [START_REF] Heise | Oxidative stress during stressful heat exposure and recovery in the North Sea eelpout Zoarces viviparous L[END_REF]. This variation was related to spontaneous, non-enzymatic GSH oxidation by emerging ROS under hyperthermia. Elevated temperatures accelerate mitochondrial respiration and increase mitochondrial ROS formation [START_REF] Heise | Production of reactive oxygen species by isolated mitochondria in the Antartic bivalve Laternula elliptica (King and Broderip) under heat stress[END_REF]. As the current summer GSH t increase was only measured at BB, it appears as an indication that glutathione pool responds to stressful conditions particularly in that area of the lagoon.

Conversely to other oxidative stress endpoints, LPO was enhanced in winter reinforcing the hypothesis that metals are inducing responses in this season.

The computed IBR values were higher in summer than winter at all sites, which is in agreement with the conclusions provided by the analysis of individual oxidative stress endpoints. Nevertheless, an identical configuration of star plots representing IBR was found in winter and summer, highlighting BB as a major impacted area in the Óbidos lagoon. Despite this clear conclusion for both seasons, IBR seemed to have a higher capacity to discern sites in winter than in summer. In fact, the estimated IBR value for BB was 100 times higher than ML in winter, whereas in summer the difference between the two sites was approximately 3 times.

Interestingly, it was in winter that spatial differences in individual oxidative stress endpoints were associated with a higher metal accumulation. In general, IBR represented a potential measure of L. aurata health at the Óbidos lagoon, providing stronger evidence of the existing impact, as previously reported in the Óbidos lagoon (Pereira et al., 2009c) and in other systems [START_REF] Beliaeff | Integrated biomarker response: a useful tool for ecological risk assessment[END_REF].

Conclusions

Results of this work provided these main findings:

-Oxidative stress responses of L. aurata gills were able to detect inter-site differences pointing BB as a critical area in the Óbidos lagoon, which was confirmed by a general stress index (IBR). Higher concentrations of accumulated Cu, Cr and Cd in winter could be on the basis of the pro-oxidant challenge at BB.

-Winter-summer variations were prevalent in gills of L. aurata, which exhibited a higher accumulation of metals in summer as well as increases of oxidative stress responses. IBR values corroborated a higher stress syndrome in summer at all sites than in winter.

Nevertheless, peroxidative damage occurred at BB only in winter despite the activation of antioxidant defences. In terms of metals accumulation, gills were more responsive in winter than summer, while oxidative stress endpoints provided evidence of existing impact at BB in both seasons.

-The applied approach combining biomarkers in gills (metal concentrations and oxidative stress) and environmental data (mainly metals and nutrients) was shown to be a useful means of assessing environmental health, particularly under a moderate contamination scenario. Gills are not only an important route of entry for contaminants, but the current study has demonstrated the ability of gills to reflect environmental status (through metal bioaccumulation and oxidative stress responses).

Figure captions Figure 1 -

 captions1 Figure captions

Figure 2 -

 2 Figure 2 -Metal concentrations in gills of L. aurata captured in winter and summer in the middle lagoon (ML) and in the inner branches (Bom-Sucesso -BS, Barrosa -BB) at Óbidos lagoon. Letters denote statistically significant differences (p<0.05) between sites (in the same survey) and between surveys (for each site): (a) vs. ML; (b) vs. BS; (w) vs. winter. (Means ± standard errors).

Figure 3 -

 3 Figure 3 -Oxidative stress responses in gills of L. aurata captured in winter and summer in the middle lagoon (ML) and in the inner branches (Bom-Sucesso -BS, Barrosa -BB) at Óbidos lagoon. Letters denote statistically significant differences (p<0.05) between sites (in the same survey) and between surveys (for each site): (a) vs. ML; (b) vs. BS; (w) vs. winter. (Means ± standard errors; mv -missing value).

Figure 4 -

 4 Figure 4 -Significant correlations between metal concentrations and oxidative stress responses in gills of L. aurata. Statistical significance and correlation coefficient are represented by p and r, respectively. Non-significant correlations are not presented.

Figure 5 -

 5 Figure 5 -Integrated biomarker response (IBR) index in gills of L. aurata captured in winter and summer in the middle lagoon (ML) and in the inner branches (Bom-Sucesso -BS, Barrosa -BB) at Óbidos lagoon.

Table 1 -

 1 Water temperature (T), salinity, dissolved oxygen (DO), ammonium (NH 4 + ), nitrate+nitrite (NO 3 -+NO 2 -), phosphate (PO 4 3-), silicate (Si(OH) 4 ) and chlorophyll a (Chl a), in winter and summer of low-tide, in middle lagoon (ML) and in inner branches (Bom-Sucesso -BS, Barrosa -BB) at Óbidos lagoon.

			T	Salinity DO	NH4 +	NO3 -+NO2 -	PO4 3-	Si(OH)4	Chl a
	Season Site	(ºC)		(%)	(µM)	(µM)	(µM)	(µM)	(µg L -1 )
		ML	9	34.6	101	2.0	9.2	0.5	6.0	1.7
	Winter	BS	11	31.2	87	18	31	2.3	20	0.7
		BB	11	28.1	96	45	40	5.6	41	1.8
		ML	20 20	34.6 34.6	95	1.0	1 1.0	2.9	5.1	4.9 4.9
	Summer	BS	23	35.2	108	1.0	1.2	2.1	4.6	3.7
		BB	24	33.3	83	0.9	1.3	6.5	12	10

Table 2 -

 2 Concentrations of metals in water at low-tide, in winter and summer, in the middle lagoon (ML) and in inner branches (Bom-Sucesso -BS, Barrosa -BB) at Óbidos lagoon. Means ± standard errors with significant differences highlighted.

			Mn		Cu			Ni	Cr		Pb	Cd
	Season Site	(µM)		(nM)			(nM)	(nM)		(nM)	(nM)
		ML	0.22±0.06	2.3±0.14	7.5±2.0	7.5±1.4	0.22±0.06	0.16±0.01
	Winter	BS	0.68±0.01	3.9±0.17	7.0±0.85 8.3±0.17 0.13±0.005	0.18±0.01
		BB	1.0±0.55		4.9±0.17	a	7.7±0.87 8.5±0.17	0.26±0.08	0.24±0.06
		ML	0.13±0.01	3.3±0.21	w	8.6±2.3	6.4±1.2	0.47±0.12 0.07±0.006	w
	Summer	BS	0.12±0.01	w	4.9±0.18	w	6.3±0.81 14±1.9	aw	0.72±0.02 0.09±0.007	w
		BB 0.04±0.02	ab	7.9±1.0	abw	8.0±0.38 14±1.3	aw	0.99±0.21	0.14±0.03

w a -significantly different from ML; b -significantly different from BS; w -significantly different from winter

Table 3 -

 3 Concentrations of metals (µg g -1 ) in suspended particles of low-tide, in winter and summer, in the middle lagoon (ML) and in inner branches (Bom-Sucesso -BS, Barrosa -BB) at Óbidos lagoon. Means ± standard errors with significant differences highlighted. -significantly different from ML; b -significantly different from BS; w -significantly different from winter

	Season Site	Mn	Zn	Cu		Ni	Cr	Pb	Cd
		ML	505±30	161±32	29±4.2		22±3.1	55±6.7	40±0.31	0.19±0.08
	Winter	BS	458±48	113±9.1	35±12		23±6.7	56±14	50±10	0.18±0.08
		BB	749±18	142±4.1 52±0.35	a	28±2.0	68±2.6	44±3.4	0.36±0.11	ab
		ML 1515±92	153±3.1	24±3.0		28±3.6	48±0.32 24±1.2	w	0.20±0.05
	Summer	BS 623±9.1 a 128±17	24±3.5		28±4.1	23±3.9	21±1.5	0.18±0.05
		BB	812±13 a 134±8.2 52±2.3 a	61±2.7 aw 51±3.3	31±1.7 b	0.18±0.01 w

a

Acknowledgements

Patrícia Pereira benefits from a PhD fellowship provided by the "Fundação para a Ciência e a Tecnologia" (FCT, Grant No SFRH/BD/17616/2004).


Angelica Ronald

Henrik Larsson

Henrik Anckarsater

Paul Lichtenstein

Running Head: A TWIN STUDY A twin study of autism symptoms in Sweden

Keywords: autism, triad, twin, behavior, genetics

This study aimed to empirically identify the number of factors underlying autism symptoms --social impairments, communication impairments, and restricted repetitive behaviors and interests --when assessed in a general population sample. It also investigated to what extent these different autism symptoms are caused by the same or different genetic and environmental influences. Autistic symptoms were assessed in a population-based twin cohort of >12,000 9-and 12-year-old children by parental interviews. Confirmatory factor analyses, principal component analyses and multivariate structural equation model-fitting were carried out. A multiple factor solution was

suggested, with nearly all analyses pointing to a 3-factor model for both males and females and at both ages. A common pathway twin model fit the data best, which showed that there were some underlying common genetic and environmental influences across the different autism dimensions, but also significant specific genetic effects on each symptom type. These results suggest that the autism triad consists of three partly independent dimensions when assessed in the general population, and that these different autism symptoms, to a considerable extent, have partly separate genetic influences. These findings may explain the large number of children who do not meet current criteria for autism but who show some autism symptoms. Molecular genetic research may benefit from taking a symptom-specific approach to finding genes associated with autism.

Autism spectrum disorders (ASDs) are currently defined by a triad of symptoms: social impairments (SIs), communication impairments (CIs), and restricted repetitive behaviors and interests (RRBIs). Family studies of individuals with ASD have reported that unaffected family members often show some autistic symptoms (termed the 'broader autism phenotype') but do not always show all three types of autistic symptoms together, that is, autistic symptoms segregate out in family members, suggesting different autistic symptoms may have different familial influences [START_REF] Bolton | A casecontrol family history study of autism[END_REF]. Recent population-based twin studies have reported that each of these three sets of features is highly heritable but appears to be caused by largely different genetic influences [START_REF] Happé | Fractionable Autism Triad': A Review of Evidence from Behavioural, Genetic, Cognitive and Neural Research[END_REF]. This evidence from family and twin studies, combined with the lack of theories in cognitive psychology that can explain all three parts of the triad together, has been cited in support of the hypothesis that autism symptoms are largely 'fractionable' [START_REF] Happé | Fractionable Autism Triad': A Review of Evidence from Behavioural, Genetic, Cognitive and Neural Research[END_REF]. It has been proposed that it might be time for researchers to give up on trying to find single explanations behind the diverse symptoms in autism and focus efforts on identifying explanations for each set of symptoms [START_REF] Happé | Time to give up on a single explanation for autism[END_REF][START_REF] Mandy | Research review: What is the association between the social-communication element of autism and repetitive interests, behaviours and activities?[END_REF][START_REF] Abrahams | Advances in autism genetics: on the threshold of a new neurobiology[END_REF].

The fractionable autism triad hypothesis is based on the idea that SIs, CIs, and RRBIs are partly independent dimensions of behavior, and ASDs occur when children show extreme forms of these problems. A contrasting hypothesis is that all three sets of autistic symptoms in the triad are part of a single underlying dimension. Many factor analytic studies have explored whether autistic symptoms fall into one or multiple statistical factors (for reviews see 2, 4). The majority of factor analysis studies report that multiple factors underlie autistic symptoms, but there are some exceptions, most notably two studies on the Social Responsiveness Scale, which report a single principal component [START_REF] Constantino | Reciprocal social behavior in children with and without pervasive developmental disorders[END_REF][START_REF] Constantino | The factor structure of autistic traits[END_REF]. A recent factor analysis of parent and teacher ratings of 730 children with Pervasive Developmental Disorder using items that directly matched the DSM criteria reported that a three-factor solution provided the best fit to the data [START_REF] Lecavalier | Deconstructing the PDD clinical phenotype: internal validity of the DSM-IV[END_REF].

To date, most twin studies on this issue have been based on a UK twin sample assessed in middle childhood at ages 7 and 8 on measures of autistic traits (with one exception, see 12). These previous twin studies have reported modest to moderate genetic overlap across different parts of the autistic triad in the general population [START_REF] Ronald | Genetic Heterogeneity Between the Three Components of the Autism Spectrum: A Twin Study[END_REF], in extremescoring groups [START_REF] Ronald | Phenotypic and genetic overlap between autistic traits at the extremes of the general population[END_REF], and in children with suspected ASDs based on a parent interview [START_REF] Dworzynski | Relationship between symptom domains in autism spectrum disorders: a population based twin study[END_REF][START_REF] Mazefsky | Genetic and environmental influences on symptom domains in twins and siblings with autism[END_REF]. These findings suggest that while some genetic influences may confer risk for developing all the symptoms of the autistic triad, others are specific to particular symptoms.

The present paper attempts to test the hypothesis of a fractionable autism triad at the phenotypic level using factor analysis, and at the etiological level using multivariate structural equation twin model-fitting. Data were collected from a Swedish population twin sample. In contrast to the previous twin studies, which used 'trait' measures, the present study employed a measure that closely reflected the current diagnostic criteria for autism [START_REF]Diagnostic and statistical manual of mental disorders[END_REF]. Our prediction, based on the results of previous factor analyses, was that we would identify three factors underlying autism symptoms. We expected to find genetic effects broadly influencing symptoms across the autism triad but also symptom-specific genetic effects.

Materials and Methods

Participants

The Child and Adolescent Twin Study in Sweden (CATSS) is a nation-wide cohort that focuses on all Swedish twins turning 9 or 12 years since July 2004. CATSS has an 80% response rate, making it a highly representative population sample. Data were available on 12446 children: N = 5944 9-year-olds and 6496 12-year-olds (the two samples were independent). Data were collected in the month of the child's birthday. 51% of the sample were male. One hundred and thirty children (76 males, 54 females) were excluded from the analyses because they had a known brain injury (N = 118) or a chromosomal syndrome (N = 12).

For 281 twin pairs, DNA was used to determine zygosity based on a panel of 48 SNPs derived for zygosity analyses. For the remaining twins, an algorithm was used based on 5 items concerning twin similarity and confusion [START_REF] Lichtenstein | The Swedish Twin Registry: a unique resource for clinical, epidemiological and genetic studies[END_REF]. In over 85% of the sample, parental interviews were conducted with the twins' mothers; no significant differences were found between mother and father ratings [START_REF] Anckarsater | Child Neurodevelopmental and Behavioural Problems are Intercorrelated and Dimensionally Distributed in the General Population[END_REF].

Measure

Parents were interviewed with the Autism -Tics, AD/HD, and other Comorbidities inventory (A-TAC; 15, 16), a telephone interview designed for large-scale epidemiological research in neuropsychiatry. Thirteen items were used to assess autistic symptoms. Items were scored "1" ("yes"), "0.5" ("yes, to some extent") and "0" ("no").

The ASD items in the A-TAC have shown good test-retest reliability (.83-.94) and have been shown to discriminate clinically diagnosed ASDs (sensitivity = .89, specificity = .78, positive predictive value = .68, area under ROC curve = 0.88) when administered by laymen over the phone [START_REF] Hansson | Psychiatric telephone interview with parents for screening of childhood autism -tics, attention-deficit hyperactivity disorder and other comorbidities (A-TAC): preliminary reliability and validity[END_REF]. All 13 items showed good internal consistency in the CATSS data (Cronbach's alpha was .81 at both ages). A total of 1.43% of the sample (78% of whom were male) were screen positive for an ASD using the A-TAC cut-off of ≥4.5 on the total score (see 17). Based on parent tick box information, a total of 0.7% of the sample (N = 87) were reported to have an ASD diagnosis and these children had elevated A-TAC scores, with mean scores of 5.52, 5.36, 4.03, and 4.44 standard deviations above the population mean for the total score, SIs, CIs and RRBIs, respectively. Together these data support the validity of the A-TAC measure and suggest that the sample is representative of the general population.

Analyses

Factor analyses. The factor structure of the thirteen ASD items was explored using principal component factor analysis. The Kaiser-Guttman criterion for factor extraction was applied, that is, number of factors was selected as the number of factors with eigenvalues of >1.

Next, confirmatory factor analyses (CFA) were conducted in the Mplus software version 4.0 to explore the adequacy of fit for the models suggested from the principal component factor analysis [START_REF] Muthén | Mplus. Statistical Analyses With Latent Variables. User's Guide[END_REF]. A robust weighted least squares estimator (WLSMV) was used, an estimation technique appropriate for categorical data [START_REF] Muthén | Mplus. Statistical Analyses With Latent Variables. User's Guide[END_REF]. Phenotypic Analyses. Phenotypic correlations were employed to explore the strength of the relationship between the autism subscales. Correlations were estimated using a saturated model in the Mx program (see below), and mean and variance differences by sex, zygosity, age and twin order were evaluated in the saturated model.

The number of children who scored at least 1 on each subscale (i.e. either at least one item was rated 'yes' or at least two items were rated 'yes, to some extent') were identified and the degree of symptom overlap in these children was explored.

The Twin Design. The twin design is based on comparing the within-pair similarity of monozygotic (MZ) and dizygotic (DZ) twins on a measure or trait of interest [START_REF] Plomin | Behavioral Genetics[END_REF]. The design is based on the assumption that MZ twins share all of their DNA and DZ twins share on average half of their DNA.

Heritability refers to the proportion of variation of a trait in a population explained by genetic influences. 'Environmental influences' in the twin design refer to all variance that is not explained by genetic influences, and is split into two types, shared and nonshared. Shared environment refers to experiences that make children growing up in the same family similar; nonshared environment refers to environmental influences that make children growing up in the same family different [START_REF] Plomin | Behavioral Genetics[END_REF].

Twin Correlations. Twin similarity coefficients (intraclass correlations; 20) were used for an initial examination of the twin data to compare MZ and DZ twin similarity, and were estimated using the Mx program [START_REF] Neale | Mx Statistical Modeling[END_REF]. A trait is under additive genetic influences (A) when twin similarity is greater for MZ than DZ twins. Shared environment (C) is indicated if DZ twin correlations are greater than half the MZ twin correlations. Nonshared environment (E) is indicated by the extent to which MZ correlations are less than unity.

When MZ twins are more than twice as similar than DZ twins, this suggests nonadditive genetic influences (D) such as dominance (interaction of two alleles at the same locus) or rater contrast effects (as described in Model-fitting Analyses).

To evaluate genetic and environmental sources of covariance across variables, cross-trait cross-twin (CTCT) correlations were calculated. CTCT correlations examine the covariance between twin 1 on variable 1 and twin 2 on variable 2, separately for MZ and DZ twins.

Model-fitting Analyses. Prior to model-fitting, scales were corrected for sex and interview order using regression. Multivariate genetic models decompose the variance of each phenotype and the covariances between phenotypes into A, C or D, and E. Three standard multivariate models were used to investigate the sources of covariance between the autism subscales: the Cholesky model, the independent pathway model, and the common pathway model. The Cholesky decomposition model is the most ubiquitous multivariate twin model; it is based on triangular decomposition and demonstrates the degree to which covariance is explained by genetic and environmental influences. The independent pathway model includes etiological influences shared between the variables as well as etiological influences specific to each phenotype. Finally, the common pathway model is the most constrained of the three and hypothesizes that the common variation between behaviors is due to a single underlying and unmeasured latent factor. This latent factor itself can be investigated in terms of genetic and environmental sources of variance.

The model also includes specific genetic and environmental sources on the individual measured phenotypes.

In accordance with the assumptions of the classic twin design, A latent variables were fixed to correlate 1.0 and .5 for MZ and DZ twins, respectively; C latent variables were fixed to correlate 1.0 for both MZ and DZ twins (because all twins in the sample were reared in the same family); D latent variables were fixed to correlate 1.0 and .25 for MZ and DZ twins, respectively; and E latent variables were fixed to correlate 0 [START_REF] Neale | Mx Statistical Modeling[END_REF].

Each of the three multivariate models was run with the following combination of variance components: ACE, ADE, AE, ADEs and AEs. "s" refers to another parameter that can be added to the model, representing a form of phenotypic interaction between twins [START_REF] Neale | Mx Statistical Modeling[END_REF]. The ACE model is a standard twin model and as such was tested, but the ACEs model was not run in addition to this because the majority of twin correlations suggested that there were no shared environmental effects and negative sibling interaction parameters were present in our data. Therefore the ACEs model was not considered a suitable model to test. When modelling parent report data, a negative phenotypic interaction often occurs because there is a contrast effect in the parental ratings of their behavior (i.e., parents inadvertently exaggerate behavioral differences between the children). Contrast effects are implicated when DZ variance is significantly greater than MZ variance.

Models were fit to raw data using the Mx structural equation modelling software [START_REF] Neale | Mx Statistical Modeling[END_REF]. This approach does not yield a χ 2 for assessing the fit of the model, however, the fit of a model can be assessed by calculating the difference between the negative log-likelihood (-2LL) of the model and that of a saturated model (i.e., a model in which the variance/covariance structure is not estimated and all variances and covariances for MZ and DZ twins are estimated). The difference in -2LL is asymptotically distributed as χ 2 with degrees of freedom equal to the difference in the number of parameters in the full model and that in the saturated model. Akaike's information Criterion (AIC; 22) and Bayesian Information Criterion (BIC; 23) fit indices were also used to compare the fit of alternative models.

Results

Factor Analysis

Principal components analysis of the thirteen A-TAC autism symptom items showed that there were consistently three factors that had eigenvalues >1. This three factor solution was found when the data were analyzed together, for males and females separately, at ages 9 and 12 separately, and with either Oblimin or Varimax rotations. The only exception was for males at age 12, for whom only two factors had eigenvalues >1.

Table 1 presents the loadings of the items on each factor using Varimax rotation on all the data (combining ages and genders). The items loading on the first factor all represented questions assessing difficulties with social interaction; this factor explained 33% of the variance in all the items. The 'make-believe' item, which falls in the CI domain in diagnostic criteria, did not load with the other communication items but loaded onto factor 1, with the four SI items. The RRBI items loaded most heavily on factor 2 (explaining 9% of variance) and the CI items loaded most heavily on factor 3 (explaining 8% of variance).

Insert Table 1 Next, confirmatory factor analysis was conducted to explore the fit of the suggested three-factor model. The fit indices suggested an excellent fit of this model (comparative fit index (CFI) = .96; Root Mean Square Error of Approximation (RMSEA) = .02). We also considered a one-factor model. The results showed that the three-factor model had a significantly better fit than this model (Δχ 2 (3, n = 12327) = 310.36.40, p < .001).

Scale Construction

Items were divided into three scales based on the results of the factor analyses. The items in the SI, RRBI and CI subscales (5 items, 5 items and 3 items, respectively) relate directly to the items in Table 1 with underlined loadings for factors 1, 2, and 3, respectively. Items in each of the scales were summed and converted into scores as a proportion of the total possible score given the number of items completed (which was required to be more than half). The internal consistencies of the SI, CI and RRBI subscales, expressed as Cronbach's alpha, were .74, .54 and .71 respectively. All scales were positively skewed and therefore inverse-transformed prior to the model-fitting. Table 2 presents descriptives for the subscales.

Insert Table 2

Phenotypic Relationship

Table 2 presents the phenotypic correlations between the three autism subscales.

All correlations were positive and significant (p<.01), ranging from .15-.42. Table 3 presents the prevalence of autism symptoms, defined as raw scores of ≥1 on each subscale, alone and in combination with other symptoms.

Insert Tables 3 and4 Data Descriptives Means could be equated without a significant decrease in the saturated model fit for twin 1 and twin 2, for males and females, and for age 9 and age 12 data, but they could not be equated across MZ and DZ twins. Variances could be equated for twin 1 and twin 2 but could not be equated across gender, zygosity or age. DZ variance was significantly larger than MZ variance for SIs and RRBIs, and MZ variance larger than DZ variance for CIs, for both genders (p<.05), and therefore sibling interaction paths were included in the models (see below).

Twin Correlations

Univariate twin correlations are presented on the diagonal for each zygosity group in Table 4. MZ twins were rated as more similar than DZ twins on all scales at both ages, suggesting significant genetic influences. In most cases, DZ correlations were half or less than half the MZ twin correlations, suggesting that there were no shared environmental influences on the subscales and that there were nonadditive genetic influences on the scales or rater contrast effects present in the data (because both have the effect of decreasing the DZ similarity to less than half the MZ similarity). MZ correlations were consistently less than unity suggesting that nonshared environment, which includes any variance due to measurement error, influenced all subscales.

Modest differences between male and female same-sex twin correlations suggested that there were some quantitative sex differences, that is, the degree to which genetic and environmental influences affect these traits may vary between the genders. Qualitative sex differences, which refer to different genetic and environmental influences affecting males and females, were not indicated by these univariate twin correlations because for the most part DZOS twin correlations were not significantly lower than DZ same-sex twin correlations.

Cross-trait cross-twin (CTCT) correlations are presented on the off-diagonal of Table 4 (male CTCT correlations below diagonal, female CTCT correlations above diagonal). Most of the MZ CTCT correlations were greater than their equivalent DZ CTCT correlations, suggesting that genetic influences were to some degree explaining the overlap between different subscales. The MZ CTCT correlations were all lower than the phenotypic correlations, suggesting that nonshared environment explained part of the covariation.

As noted above, the twin correlations did not suggest that qualitative sex differences were present in the data, and significantly different MZ and DZ variances suggested that sibling interaction paths were required in the model. Therefore, in line with previous behavior genetic research, DZOS were excluded from the models because of the presence of sibling interaction paths (it is considered too complex to include both qualitative sex differences and multiple sibling interaction paths in the same multivariate model, see [START_REF] Simonoff | Genetic influences on childhood hyperactivity: Contrast effects imply parental rating bias, not sibling interaction[END_REF].

Insert Table 5

Model Comparisons

Table 5 presents the fit statistics for models at both ages. In all these models, means were equated across twins 1/2 and males/females but not for MZ/DZ groups. Fit statistics for the ACE, ADE, AE versions are available from the first author on request

The best-fitting model was established based on three fit indices: the change in the Chi-square between submodels (smaller changes with the highest number of degrees of freedom were favored), the lowest AIC, and the most negative BIC value. First, ACE and ADE Cholesky models were run in which all parameters were estimated. Then the AE model was tested (by dropping the C/D paths) and there was no significant decrease in fit compared to the ACE model (Δχ2 = 6.34, (12 df), p=n.s.) or the ADE model (Δχ2 = 17.58, (12 df), p=n.s.). Sex differences could not be dropped from this model without a significant decrease in fit (Δχ 2 = 216.08 (12 df), p<.001). Next, independent and common pathway AE models with sex differences were compared to the Cholesky model. The AE common pathway gave the best fit compared to the Cholesky and independent pathway models.

Models were run with sibling interaction paths and there was a deterioration in fit when we tested a model without sibling interaction paths (Δχ 2 = 39.05 (3 df), p <.001 ), but the sibling interaction paths could be equated across gender without a significant deterioration in fit (Δχ 2 = 1.03 (3 df), p=n.s.).

For the 12-year-olds, very similar results were found and the best-fitting model was also the common pathway AEs model with different estimates for male and females and sibling interaction paths that could be equated across gender.

Testing the Significance of Symptom-Specific Genetic Effects

The lower section of Table 5 shows the fit statistics for the common pathway AEs at both ages when the paths for specific additive genetic influences on each subscale were sequentially dropped from the model. It was not possible to drop the paths for specific additive genetic influences for SIs, CIs, or RRBIs without a significant decrease in fit in the age 9 data (SIs: Δχ 2 = 34.93 (2 df), p < .001; CIs: Δχ 2 = 82.51 (2df), p = <.001; RRBIs: The heritability of each scale can be derived from summing the percent variance explained by genetic effects in common with other symptoms with the percent variance explained by unique genetic effects in Figure 1. SIs in 9-year-old boys, for example, have a heritability of 31%+36% = 67%. Heritabilities for all symptom scales in the 9 and 12year-old children ranged from 49% (CIs in 9-year-old girls) to 76% (SIs in 12-year-old boys). The remaining variance for all scales was accounted for by nonshared environment.

Δχ 2 = 27.
In 9-year-olds, sibling interaction paths were estimated at -.14, .05 and -.04 for SIs, CIs, and RRBI subscales, respectively, and in 12-year-olds as -.15, -.02 and -.09 for SIs, CIs, and RRBIs, respectively.

In terms of the covariance between subscales, the models in Figure 1 1 these values are presented as a percentage of the total variance. Thus, for SIs for 9-year-olds, heritability was 67%, and 47% of this genetic variance was specific. For CIs and RRBIs an even larger proportion of the genetic influences were specific to that subscale. For all three subscales, the larger part of the nonshared environmental effect was specific to each subscale (66-96%). Similar results were found at age 12, as shown in Figure 1.

Discussion

This study explored the factor structure of autism symptoms in the general population in 9-and 12-year-old children. At both ages, ASD symptoms split into three factors -social impairments (SIs), communication impairments (CIs) and restricted repetitive behaviors and interests (RRBIs) -, which mirrored the theoretical autism triad as well as previous empirical results of autism symptoms in children with pervasive developmental disorders [START_REF] Lecavalier | Deconstructing the PDD clinical phenotype: internal validity of the DSM-IV[END_REF]. The only divergence from the DSM-IV division of items was the 'make-believe' item, which loaded with SIs rather than CIs (the same was also found elsewhere, [START_REF] Lecavalier | Deconstructing the PDD clinical phenotype: internal validity of the DSM-IV[END_REF]. Our results showed that many more children show one part of the ASD triad than all three together.

The twin analyses revealed that SIs, CIs and RRBIs are all highly heritable (49-72%), which agrees with most previous studies of middle childhood onwards, [START_REF] Ronald | Genetic Heterogeneity Between the Three Components of the Autism Spectrum: A Twin Study[END_REF][START_REF] Hoekstra | Heritability of autistic traits in the general population[END_REF][START_REF] Constantino | Intergenerational transmission of subthreshold autistic traits in the general population[END_REF][START_REF] Constantino | Autistic traits in the general population: a twin study[END_REF][START_REF] Constantino | Genetic structure of reciprocal social behavior[END_REF][START_REF] Ronald | A twin study investigating the genetic and environmental aetiologies of parent, teacher and child ratings of autistic-like traits and their overlap[END_REF][START_REF] Ronald | The Genetic Relationship Between Individual Differences In Social And Nonsocial Behaviours Characteristic Of Autism[END_REF][START_REF] Skuse | Measuring autistic traits: heritability, reliability and validity of the Social and Communication Disorders Checklist[END_REF][START_REF] Edelson | Genetic and environmental influences on autistic-like behaviors in 2-year-old twins[END_REF]; a lower heritability estimate has been reported in one twin study of 2-year-olds [START_REF] Constantino | Genetic structure of reciprocal social behavior[END_REF]. The ASD subscales shared a substantial degree of genetic influences, but it was also found that each ASD subscale had significant genetic influences that were specific to itself.

Similar to other studies, environmental effects were primarily nonshared and subscale-specific. With a few exceptions [START_REF] Ronald | A twin study investigating the genetic and environmental aetiologies of parent, teacher and child ratings of autistic-like traits and their overlap[END_REF][START_REF] Edelson | Genetic and environmental influences on autistic-like behaviors in 2-year-old twins[END_REF], these results concurred with previous twin studies in finding negligible shared environmental effects. Small differences were found between ages 9 and 12, for example, heritabilities were on average slightly higher at age 12 than at age 9, but overall the results are noticeable for their similarity across ages.

The best-fitting model included different parameter estimates for males and females, in line with previous studies (e.g. ( 9), but overall sex differences were modest.

Previous results [START_REF] Ronald | Genetic Heterogeneity Between the Three Components of the Autism Spectrum: A Twin Study[END_REF][START_REF] Ronald | Phenotypic and genetic overlap between autistic traits at the extremes of the general population[END_REF][START_REF] Dworzynski | Relationship between symptom domains in autism spectrum disorders: a population based twin study[END_REF][START_REF] Ronald | The Genetic Relationship Between Individual Differences In Social And Nonsocial Behaviours Characteristic Of Autism[END_REF] and the present study concur in suggesting that there are a degree of overlapping genetic and environmental influences between different autism behaviors in the general population but also significant genetic and environmental influences that are specific to each part of the autism triad. One slight difference is that while the results from the UK twin studies found that SIs and RRBIs showed the least amount of phenotypic and genetic overlap within the triad, the present study found that this overlap was highest between SIs and RRBIs. The methodology and sampling were similar across the two studies. This difference therefore may be a result of the different measures employed: the A-TAC is a shorter measure that aims to directly reflect ASD symptoms, whereas the CAST, used in the previous UK twin studies, has a greater number of items that may be considered more 'trait-like'.

The present study should be considered in light of its limitations. Although the items used here closely mirrored the DSM-IV criteria, parental interviews are not equivalent to a psychiatrist's diagnosis (which would not be feasible with a population cohort) and the A-TAC is not currently considered a gold standard questionnaire for diagnosis in the same way as, for example, the longer and more established ADI-R. The scales had skewed distributions, which may lead to bias in parameter estimates if the assumption of multivariate normality is not met [START_REF] Derks | Effects of censoring on parameter estimates and power in genetic modelling[END_REF]. Possibly because of the skewed distributions, all the models fit significantly worse than the saturated model, which is another limitation of the data. Parent report contains some bias [START_REF] Najman | Bias influencing maternal reports of child behaviour and emotional state[END_REF] and shows only modest correlation with other raters when studying autistic symptoms [START_REF] Ronald | A twin study investigating the genetic and environmental aetiologies of parent, teacher and child ratings of autistic-like traits and their overlap[END_REF]. Parental assessment of problem behavior is a practical option for large studies, and parents are familiar with behavior across time and a range of situations. Finally, diagnostic criteria are likely to evolve, and it will be important to explore the genetic architecture of autism behaviors in more detail, for example, including RRBI subcategories and additional symptoms such as hypersensitivity.

This study has important implications for clinicians. First, this study reports that many children showed autism symptoms in part of the triad only, in as severe a form as might be expected in children with a formal diagnosis of autism or Asperger's disorder.

Further research needs to follow up these children in more detail, to explore their clinical needs and underlying cognitive deficits; some may be 'phenocopies' in that their problems are due to different underlying disorders such as anxiety. Further research could address whether RRBIs on their own are more common that SIs or CIs on their own, as suggested by the results in Table 3, or if these different frequencies are an artifact of the measuring instrument. Nevertheless, this finding appears to concur with the high prevalence of Pervasive Developmental Disorder-Not Otherwise Specified (PDD-NOS) diagnoses, which make up almost a third of all ASD diagnoses and are given when presentations do not meet the criteria for Autistic Disorder because of late age at onset, atypical symptomatology, or subthreshold symptomatology, or all of these. Diagnostic criteria that attended to the level of impairment of a child on each aspect of the autism triad, rather than just considering a total sum of autistic symptoms, might better classify these heterogeneous clinical groups.

A child with ASD with predominantly social impairments could have different educational and clinical needs as compared to a child with predominantly communication and repetitive behavior problems.

Molecular genetic research has begun to explore the possibility of symptomspecific genetic influences in autism, using candidate gene approaches [START_REF] Brune | 5-HTTLPR Genotype-Specific Phenotype in Children and Adolescents With Autism[END_REF], linkage [START_REF] Liu | Genome-wide linkage analyses of quantitative and categorical autism subphenotypes[END_REF] and genome-wide association [START_REF] Ronald | genome-wide association study of social and non-social autistic-like traits in the general population using pooled DNA, 500K SNP microarrays and both community and diagnosed autism replication samples[END_REF]. Complexity continues to be a key feature in the molecular genetics of ASD; the present findings suggest that knowledge about the causal pathways underlying individual symptoms will help to lead the way in this area of research.

17

 17

1 Figure 1

 11 Figure 1 presents the unsquared unstandardized path estimates for the best-fitting

  divide genetic and environmental influences into those that are common to the three subscales, shown in the top half of each model influencing the common latent factor (L), and those that are specific to each subscale individually, shown in the bottom part of each model. For the 9-year-olds, using the path diagram, the percentage of genetic influences for SIs that were specific to SIs in males (i.e., not in common with the genetic influences on CIs and RRBIs) was calculated as follows: ((.59*.59)/ ((.59*.59)+(.72*.87*.87*.72)))*100% = 47%. Conversely, the remaining genetic influences i.e. ((.72*.87*.87*.72)/ ((.59*.59)+(.72*.87*.87*.72)))*100% = 53% of the genetic influences on SIs, were shared with CIs and RRBIs via the common latent factor (L). In the bar chart in Figure

  Note. Sex = Quantitative sex differences. Sibling interaction paths (s) were equated across gender. a Best-fitting model at both ages. -2LL = log likelihood fit; df = degrees of freedom; par = parameters; LRT(df) = likelihood ratio χ 2 test with Δ df comparing model to either saturated model (Model 1) or to best-fitting model (Model 4); AIC/BIC = Akaike's/ Bayesian Information Criterion. A/D = Additive/ nonadditive genetic influences, E = Nonshared environment.

Figure 1 .

 1 Figure 1. Path diagrams of best-fitting common pathway model in 9-year-olds (top panel)

Table 1 .

 1 Loadings of autism symptom items in principal component analyses.

	Factor/Item

Table 2 .

 2 Descriptives and phenotypic correlations Note. SIs = social impairments; CIs = communication impairments; RRBIs = restricted repetitive behaviors and interests. Male phenotypic correlations presented below diagonal, female above diagonal. All phenotypic correlations p<.01.

	SIs	CIs	RRBIs
	M	0.11	0.11	0.24
	SD	0.43	0.35	0.59
	Skewness	6.06	3.95	3.65
	Range	0-5.00	0-3.00	0-5.00
	Age 9 (N = 1517/1426)
	SIs	-	.21	.33
	CIs	.27	-	.14
	RRBIs	.37	.23	-
	Age 12 (N = 1651/1566)
	SIs	-	.30	.40
	CIs	.31	-	.24
	RRBIs	.41	.27	-

Table 3 .

 3 Prevalence of autism symptoms alone and in combination with other symptoms in the sample. Note. Sample consisted of one twin per pair selected randomly for birth order. Presence of autism symptom defined as a raw score of 1 or above on a subscale (see text). SIs = Social impairments; CIs = Communication impairments; RRBIs = Restricted repetitive behaviors and interests.

		Males (N = 3165) Females (N = 2991)
	SIs only	41 (1.3%)	48 (1.6%)
	CIs only	155 (4.9%)	79 (2.6%)
	RRBIs only	284 (9.0%)	166 (5.5%)
	SIs + CIs, No RRBIs	17 (0.5%)	14 (0.5%)
	SIs + RRBIs, No CIs	62 (2.0%)	40 (1.3%)
	CIs + RRBIs, No SIs	56 (1.8%)	18 (0.6%)
	SIs + CIs+ RRBIs	63 (2.0%)	11 (0.4%)

Acknowledgements

The authors would like to thank the participants of CATSS. This study was supported in part by the Swedish Council for Working Life and Social Research and the Swedish Research Council. HL was supported by a postdoctoral stipend from the Swedish Brain Foundation and the Karolinska Institutet Center of Neurodevelopmental Disorders, Stockholm. AR was funded by the Royal Society. The authors report no conflicts of interest.