**Note that below we have provided additional information for our BSI proposal to aid reviewers. We provide a table of contents for you to quickly view a particular section or you may simply scroll through this page.

www.utexas.edu/tmm/tnhc/fish/research/truchas_mexicanas/

or

www.americanfishes.com/mexico

Click below to view

Accomplishments from PI Mayden's Tree of Life Grant (ATOL), years 1 and 2

Table of Contents

Collaborators of Truchas Mexicanas
Terrain of the Sierra Madre Occidental (SMO) - Habitat of Truchas or Aparaques
Devestating Impacts on the Aquatic Ecosystems of the SMO
Sampling and Specimen Processing in SMO
Ecological Niche Modeling -- GARP Modeling Success
New GARP Methods
Morphological Analyses and Taxonomic Study of Mexican Trout
     Shape differences
     Meristic analyses
     Quantitative analyses of spotting and pigmentation differences
Molecular Analyses for Population Genetics and Phylogeny of Mexican Trout
     Genetic Studies
     Observed Microsatellite Variation: Preliminary Data
     Observed DNA Sequence Variation: Preliminary Data
Phylogeny & Zoogeography
Comparative Phylogeography
Amphibian Biodiversity in Sierra Madre Occidental
eLife Electronic Data System for Research and Education
Literature Cited

Collaborators of Truchas Mexicanas over several years

Current Collaborators of Truchas Mexicanas

Return to Table of Contents

Terrain of the Sierra Madre Occidental (SMO) - Habitat of Truchas or Aparaques

The following photographs depict some of the general terrain where we are conducting our research efforts for the aquatic fauna of SMO. It is quite rugged and requires considerable time for most access. More photos are available at http://www.americanfishes.com/mexico/index%20photo%20gallery.htm

Return to Table of Contents

Devestating Impacts on the Aquatic Ecosystems of the SMO

As outlined in the proposal there are many critically important environmental impacts on the aquatic and terrestrial ecosystems of the SMO through human activities and the introductions of cultured Rainbow Trout. Below, we provide a few images to highlight these activities.

Typical aquaculture facility literally littering the SMO and producing Rainbow Trout

Left, Rainbow Trout growout facility constructed by Mexico; Right, make-shift growout facility produced by local resident. The latter is a diversion of a local stream and contains Rainbow Trout that will be ultimately released as feral trout in the stream after major rains and high water.

Livestock grazing (horses and cattle) that denude the landscape and lead to high levels of erosion, especially near stream banks where they congregate for water.

Left, Logging of major trees along margins of rivers and streams and throughout watershed result in significant erosion into stream habitats. Right, farming of near-stream areas with no erosional barriers result in high sediment loads into the streams during the rainy season and can eliminate native habitats for fishes, amphibians, and macroinvertebrates (obvious food source for the former groups)

Return to Table of Contents

Sampling and Specimen Processing in SMO

At each location USA and Mexican team members collect samples of fishes and often macroinvertebrates and amphibians, as well as extensive locality, environmental, and field data.

Return to Table of Contents

Ecological Niche Modeling -- GARP Modeling Success!

GARP modeling is essential for exploration of difficult terrain; we have demonstrated the success of this method in both laboratory analyses and in the field!! This is detailed below.

FIGURE. Map of rivers of northern Mexico illustrating our successful use of GARP in predictions of occurrences of native trout in Mexico. Green rivers = Model convergences and high-probability predictions of appropriate "niche" habitats specific for trout in Mexico. Red rivers = low probability of trout habitat. Black dots represent data used to develop model predictions depicted by green rivers. Note that the predicted green rivers overlay the ranges of the different trout species we have discovered with field surveys, including the upper Conchos River system (Atlantic Slope drainge) and southern Pacific drainagse where trout have never been discovered before.

FIGURE. Enlarged map of rivers of northern Mexico illustrating our successful use of GARP in predictions of occurrences of native trout in Mexico. River colors as above. Inset legend explains taxonomic forms of trout in this study and their drainages. High probability predictions of habitat (green rivers) based only on Rio Yaqui samples (black dots) provide excellent and highly significant (P < 0000.1) predictions of native trout in Mexico. Note green rivers in the Rio Conchos basin (screened pink) and further south in the headwaters of the Rios Nazas and Mezquital, providing highly significant predictions that native trout exist in these basins but have never been discovered. The first occurrence of ANY native rainbow lineage trout in an Atlantic basin drainage is in the Rio Conchos (orange trianges), predicted to occur there by GARP analysis and discovered in 2005.

Georeferenced localities of trout used to examine distributions and develop niche models are extracted from museum specimen records and personal collections, scientific literature and our field surveys (Hendrickson et al. 2003 and additional unpublished data).  Environmental data sets in our initial modeling included 13 digital maps ("coverages") summarizing aspects of topography (elevation, topographic index, slope, aspect, and flow accumulation, from the US Geological Survey's Hydro-1K data set, native resolution 1x 1 km), climate (as absolute and average minimum and maximum temperature (1:1,000,000), precipitation (1:1,000,000) and evapotranspiration (1:4,000,000) polygon coverages) and soils (1:250,000 polygon coverage) (climate and soil coverages from the CONABIO digital map library, available at: http://conabioweb.conabio.gob.mx/metacarto/metadatos.pl).  All coverage data were resampled at 0.01 x 0.01 decimal degree resolution and clipped into 2 datasets; one consisting of all of México and one to a buffered distance of 100 km of the locality data.  Models were built on the buffered set and then projected onto all of México to facilitate analysis. 

GARP (Genetic Algorithm for Rule-set Prediction; Stockwell and Peters, 1999; Anderson et al., 2003; Peterson and Robins 2003) is a previously tested approach for estimating species' "ecological niches" or their distributions. GARP is a machine-learning approach for detecting associations between known occurrences of species and landscape characteristics, particularly when the ecological landscape is complex and relationships are not straightforward (Stockwell and Noble, 1992; Stockwell and Peters, 1999).  Inspired by models of genetic evolution, GARP models are composed of sets of rules that "evolve" through an iterative process of rule selection, evaluation, testing and incorporation or rejection (Holland, 1975; Anderson et al., 2002). GARP searches for non-random associations between occurrences of species (georeferenced localities in geographic coordinates) and environmental variables (i.e., digital maps of relevant ecological parameters).  GARP divides occurrence data randomly into data input into the genetic algorithm for model development (“training”) and an independent data set ("extrinsic testing data") for evaluation of model quality at user-determined proportions (80% and 20%, respectively). The training data are resampled to 1250 training data points, and 1250 points are randomly sampled from all cells where training data do not occur (as a pseudo-absence data set).  These 2500 points are used to develop the rules and as intrinsic testing data (for evaluation of rule predictivity).  Spatial predictions of presence and absence can hold two types of error: omission (areas of known presence predicted absent) and commission (areas of known absence predicted present), which can be summarized in a measure of predictive accuracy as the percentage of points correctly predicted as present or absent (the correct classification rate of Fielding and Bell (1997)). Changes in predictive accuracy from one iteration to the next are used to evaluate whether particular rules should be incorporated into the model or not, and the algorithm runs either 1000 iterations or until convergence (Stockwell and Peters, 1999). The final ecological niche model rule-set is then projected onto the digital maps that are the environmental data sets input into the algorithm to identify areas fitting the model parameters, a hypothesis of the potential geographic distribution of the species. Predictive accuracy of the resulting model is evaluated based on the extrinsic testing data.  Accuracy metrics based on these data are used to select the 10 best-subset models from all models generated.

New GARP Methods

GARP methods have evolved since we modeled and successfully located trout populations not previously known in Mexico. In our future analyses we will modify the methods employed above to include the following to more accurately predict the occurrences of native trout in Mexico

The environmental data sets to be included in future analyses will be 45 digital maps (”coverages”) summarizing variables used in earlier analyses and other variables (e.g., flow direction and flow accumulation, 16-day composite remotely-sensed data layers (one composite per month during 2003 of the Normalized Difference Vegetation Index (NDVI) and the Enhanced Vegetation Index (EVI), from the NASA-MODIS/Terra data set;  native resolution 500 x 500 m; Justice et al., 1998). Vegetation indices are numerical quantities derived from spectral reflectance values that are intended to be proportional to a certain geo-biophysical variable. In our case, both indices are derived form reflectances measured in the visible and near-infrared domains and they are complementary--while NDVI is sensitive to photosynthetic activity (Tucker, 1979), EVI (L = 1; C1 = 6; and C2 = 7.5 as equation coefficients) is more responsive to canopy structural variations (Running et al., 1994; Huete et al., 1999, 2002). All geographic data were resampled at 1 km resolution to facilitate analysis.

Our new GARP analyses will divide occurrence data randomly into data input into the genetic algorithm for model development and an independent data set (“extrinsic testing data”) for evaluation of model quality at user-determined proportions (50% and 50%, respectively). The input data are further subdivided into training data (for rule development), and intrinsic testing data (for evaluation of rule predictivity). Predictive accuracy of rules are evaluated based on the extrinsic testing data. Spatial predictions of presence and absence can hold two types of error: omission (areas of known presence predicted absent) and commission (areas of known absence predicted present), which can be summarized in a measure of predictive accuracy as the percentage of points correctly predicted as present or absent (the correct classification rate of Fielding & Bell (1997)). Changes in predictive accuracy from one iteration to the next are used to evaluate whether particular rules should be incorporated into the model or not, and the algorithm runs either 1000 iterations or until convergence (Stockwell & Peters, 1999). The final ecological niche model rule-set is then projected onto the digital maps that are the environmental data sets input into the algorithm to identify areas fitting the model parameters, a hypothesis of the potential geographic distribution of the species.

Given the stochastic nature of GARP (both via sampling of occurrence data and the genetic algorithm itself), GARP produces distinct results for alternate runs of the same input data, representing alternative solutions to the computational challenge. Following recently proposed best-practices approaches (Anderson et al., 2003), we will develop 100 replicates of each model; of these, we will retain the 20 models with lowest omission error, and then discard the 10 models of these 20 that present the most extreme values of area predicted present (the commission error index of Anderson et al. 2003). These “best subset” models were summed to produce final predictions of potential distributions.

Models were validated using Receiver Operating Characteristic (ROC) analysis, which evaluates model performance independently of arbitrary thresholds at which presence might be accepted (Zweig and Campbell, 1993; Fielding and Bell, 1997; Manel et al., 2001; Wiley et al., 2004). ROC assesses model performance by plotting sensitivity (proportion of presence points correctly predicted) vs. 1 - specificity (proportion of absences correctly predicted) across all possible thresholds. Because true absences are not available for all areas, pseudoabsences (any pixel without confirmed presence) are employed as surrogates (Wiley et al., 2004). Good model performance is characterized by large areas under this curve (AUC; maximizing sensitivity for low values of 1 - specificity), and z-values above the critical level. The AUC can be translated into the probability that a model is predicting presence better than at random: AUC values of 0.5 indicate no accuracy (equivalent to the performance of a random model), values 0.5-0.7 indicate low accuracy, values of 0.7-0.9 indicate useful applications, and values of >0.9 indicate high accuracy (Hanley and McNeil, 1982; Manel et al., 2001). A z-statistic can be used to compare observed AUC with the random AUC, or between AUCs for two independent analyses; if z > 1.96, then the probability is < 0.05 that the observed difference would be expected at random.

The overall approach used to compare ecological niches of trout in different regions will consist of three steps. (1) We will use a modification of the jackknife approach of Peterson and Cohoon (1999) to determine the most suitable combination of environmental data sets for predicting the species’ geographic distribution. (2) We will test the accuracy of predictions of geographic distributions using occurrence data from the region upon which the model was based. Finally, (3) we will project the ecological niche model for each area onto each other area and tested predictive ability to assess the strength of niche differences.

The first step is designed to optimize the environmental data set for a particular situation. Hence, we initially use the range of a particular species and the Mexico occurrence data to develop models. Ecological niche models based on occurrences from the native ranges are used to predict Mexican ranges, and vice versa. Multiple predictions will be developed for each occurrence data set, representing all combinations of suites of environmental data sets. Based on results of these preliminary analyses (using the ROC procedure described above), we chose suites of environmental data sets for further analyses.  For testing (step 2), when possible given sample sizes available, for validations within regions, occurrence data are divided into two categories for predicting and testing model performance, as follows. Native-range occurrences of the different species are divided into 2 x 2º checkerboards. We validate within-region model predictivity by using one of the categories listed above to predict the other category, and vice versa, and evaluated model performance using ROC. Finally, (3) each regional occurrence data set within and between species are modeled using GARP, and the resulting ecological niche model projected to each other area. Model performance in predicting ‘other’ regions is evaluated using the ROC procedure described above.

Because the terrain of the Sierra Madre Occidental is very rugged and time consuming to navigate, we sought to optimize our collecting efforts of communities likely to harbor native trout; GARP analysis proved to be extremely beneficial in not only saving valuable time but in also discovering new diversity from rivers not previously known to contain this biodiversity.  Furthermore, because we had evidence from historic correspondence (reviewed in Hendrickson et al. 2003) of native trout in the Río Conchos basin, and our review of drainage patterns, faunal distributions and habitats in the upper Río Conchos indicated that trout could very well exist in these headwaters, we sought the use of GARP analysis to predict most probable streams in this system for maintaining trout populations.  For these analyses we used environmental parameters derived from remotely-sensed imagery, interpolated station data, and and known trout habitats in the adjacent Río Yaqui and other basins. GARP determines characteristics of a species native habitat based on environmental parameters observed at known collection localities. Similarities across collection locations are used to form conclusions about the range of the species based on factors such as temperature, precipitation, elevation, soil type, etc. As there is a stochastic element to the GARP process, there is no unique solution.  Consequently, a number of models are generated (300 in this case) and a best-subset is chosen. By intersecting the 10 best-subset model convergence map with a stream coverage, streams could be assigned respective model convergence values. These values can then be ranked to reflect likelihood of native trout populations — based upon information known about known localities of native trout in the Sierra Madre Occidental. Using this algorithm we generated a map of likely locations for trout throughout the Sierra Madre Occidental, including the Río Conchos. The model correctly identified areas known to contain trout both in the remaining parts of Pacific-slope drainages, and predicted trout occurrences in tributaries of the upper Río Conchos.

A prediction of the distribution of Trout in the SMO is illustrated below in a 10-model convergence intersected with streams in Mexico. Black dots are real localities for the Yaqu Trout, some of which were used in model construction; others were used to test models. Red rivers indicate low model convergence (or probability of trout "niche" habitat"), orange is medium prediction, and green streams are those with high probability of containing trout.

Return to Table of Contents

Morphological Analyses and Taxonomic Study of Mexican Trout

Trout Species in Mexico. To date, our combined morphological and molecular investigations have identified up to ten different forms of trout in Mexico that we currently believe to be distinct species. This diversity is based on both diagnostic character arguments as well as phylogenetic placement of these taxa relative to the Mexican Golden Trout (O. chrysogaster), Nelson's Rainbow Trout (O. nelsoni; Baja California), and the Apache (O. apache) and Gila trouts (O. gilae). A separate webpage is linked here for reviewers to see a preview of the diversity and some of the morphological (coloration) diagnostic traits)

Shape differences: Geometric morphometric analyses will be used to quantify shape differences between populations of Oncorhynchus. Images of lateral and dorsal aspects of preserved specimens are against a photo-grey background and ruler. We employ a Truss with 20 mid-sagittal homologous landmarks, some from our previous analysis (Ruiz-Campos et al. 2003) (see Figure below). Cartesian coordinates of relative positions of landmarks are captured using the tpsDig (Rohlf, F.J., unpubl.). All analyses are performed separately for sexes. Centroid size are used as an estimator of body size (sensu Bookstein 1991), while body shape will be quantified using shape coordinates. Shape variables are derived from landmark data using a generalized Procrustes analysis. Two landmarks are chosen as baseline landmarks, and the relative position of the remaining landmark constellation transformed into shape coordinates (Bookstein 1991). Allometric variation is assessed both by pairwise linear regressions of ln centroid size onto ln of the character of interest, as well as by multivariate regression of shape onto centroid size (using tpsRegr; Rohlf, F.J., unpubl.). With significant allometry samples are partitioned into size classes for analyses. Size variation is evaluated with ANOVA to test for differences in centroid size among populations. Shape differences are compared visually by plotting means of shape coordinates for each population, statistically assessed with MANOVA for each of the non-baseline coordinate pairs, and summarized using CVA (Marcus 1989). Classification from CVA of shape coordinates are compared to those of 100 randomized permutations of the data. We will use relative warp analysis (RWA) to describe major trends in shape variation (Bookstein 1991; Rohlf 1993). RWA are principal components of the degree of deformation needed to transform landmarks from each specimen to a tangent configuration (Rohlf 1993). We will compare overall patterns of morphological shape variation by deriving pairwise distances among all populations from means for relative warp axes. Visual inspection of warp grids provide an intuitive means of assessing differences in shape across samples.

Figure. Lateral view of truss and standard measurements used for species of trout to evaluate shape. Truss and standard measurements for doral and ventral aspects not shown.

Meristic analyses: Considerable information on meristic variation in Oncorhynchus is available. Standard counts are utilized, to allow comparison with previously published works. These include counts of scales, gill rakers, fin rays and vertebrae (x-rays). The Wilcoxon two-sample test will be used to identify differences between sexes or OTUs for meristic characters. Meristic variation between OTUs is summarized using PCA (correlation matrix). If sexes are significantly different they are analyzed separately. Since each ANOVA will include more than two populations, the alpha value is adjusted with a sequential Bonferroni correction (Hochberg 1988) to account for multiple comparisons. As this approach may artificially inflate the alpha value, the Tukey HSD test will also be used and the most conservative estimate employed.

Quantitative analyses of spotting and pigmentation characteristics.--Prior studies have noted considerable variation in spotting patterns. Characteristics to be examined include total number of spots on dorsum, relative distribution of spots (number of spots in each of four quadrants between the lateral line and the dorsal midline, and four quadrants between the lateral line and the ventral midline), mean size of spots in each quadrant, number and area of primary and secondary parr marks, number of spots in each lobe of the caudal fin and in the dorsal fin. ImageTool v3.0 (IT; C.D. Wilcox, S.B. Dove, W.D. McDavid and D.B. Greer; The University of Texas Health Science Center, San Antonio, Texas) will be used to extract quantitative data about colors, spots (size, number, etc.) from images of both live and preserved specimens. Some color is lost upon preservation, while others are more apparent after contact with fixative (melanophores that are covered by silvery guanine). IT supports scripting and allows automated image analysis, as well as spatial calibration for linear and area measures. Density calibration will be done relative to radiation (for color) or optical density (for melanophores) standards.

Return to Table of Contents

Molecular Analyses for Population Genetics and Phylogeny of Mexican Trout

Genetic Studies—Genetic studies include microsatellites (microsats) and DNA sequencing (sequencing). We will continue to examine 12 informative microsat loci loci (Omy2, Omy27, Omy77, Omy207, Omy325, Onem8, Onem11, Ots1, Sfo8, Ssa14, Ssa85 and Ssa289) to investigate species boundaries, comparisons of variability of hatchery and native trouts, and investigate areas of introgressive hybridization. Microsat loci and methods of analyses are outlined by Nielsen and Sage (2001) and members of our team have not only trained in Nielsen’s lab for these but have continued these exact analyses in the last 4 years in the Mayden and Garcia de Leon (see letters) labs. We have budgeted to examine 1500 individuals for these loci; these data will be collected by Mexican collaborators Dr. Francisco Garcia de Leon and Anabelia De los Santos Camarillo in La Paz, Mexico; graduate student Anabelia will be working in Mayden’s lab Aug.-Dec. 2006 on data generation and analysis of some hatchery and some native trout that they and Dr. Barriga have sampled from different locations in Mexico; these data will serve as critical comparative allelic information for comparisons with native trout in hybrid analysis (PI Mayden is providing financial support and laboratory space and training to accommodate this analysis and Dr. Irene Barriga Sosa from Univ. Autó. Metro., México, is providing support in this collaborative effort from her SAGARPA-CONACYT-2005-12147 award to examine genetics of hatchery trout, see letters). In addition to Mayden’s lab working with microsats of these tetraploid species, Mayden is also in collaboration with Dr. R.M. Wood (also at SLU) in microsat analyses of tetraploid sturgeon of Scaphirhynchus; thus PI Mayden has demonstrated both laboratory and analytical experience with various microsat applications in tetraploid fishes.

Total genomic DNA will be extracted with QIAGEN DNEasy kits following QIAGEN protocols. Twelve microsatellite loci will be amplified as in our preliminary studies using standard PCR protocols described in Nielsen and Sage (2001); all microsatellite data generation for prelimary studies was conducted in Dr. Jennifer Nielsen’s laboratory, Anchorage, Alaska, to ensure that are data are comparable to her previously published data. Microsatellites will be visualized on a LI-COR Long Reader 4200 automatic sequencer with Gene ImagIR (v3.0) software. Mitochondrial D-loop sequence data will be generated amplified using primers DL1F and DL4R (Bagley and Gall 1998). ATPase 6 and 8, and GH1C and GH2C will be amplified using primers from Kontula et al. (2003) and Phillips et al. (2004). QIAGEN gel purification kits will be used to purify PCR products prior to sequencing on an ABI3100 Genetic Analyzer. Complete sequences will be assembled from individual electropherograms, edited and aligned with BioEdit v6.0.7 (Hall 2001).

Tests of linkage disequilibrium, HWE and genotypic differentiation will be performed using GENEPOP (Raymond and Russet 1997). Genetix will be used for Factorial Correspondence Analysis. Analysis of molecular variation (AMOVA) will be performed using Arlequin v2.0 (Schneider et al. 2000).

We will also generate sequences of mtDNA and nDNA genes (Control Region [CR], ND2, GH, ATPase 6/8 = 5090 bp of mtDNA data per specimen) for 150 individuals per year. To optimally identify individuals we will use microsat data as a rough screen of lineage diversity, and choose individuals representing all lineages. CR has already been surveyed for much of the range of rainbow trout (Bagley and Gall 1998) and a small subset of populations in the northern SMO and Mexican Golden Trout (Nielsen et al. 2001) and found to be highly informative for delineation of taxa and further has extremely interesting indels that may have functional significance in the adaptation of these southern trout (Coskun et al. 2003). These genes have also been used in several studies by Dr. Ruth Phillips and collaborators (Domanico and Phillips 1995, Domanico et al 1997, Phillips and Oakley 1997, Phillips et al 1995, Reed et al 1998a,b). We have continued to examine the CR for both aspects and have found that it provides phylogenetic resolution for materials we have (Figure below in preliminary data). Nuclear genes include well-characterized introns in a duplicated growth hormone, GH1 and GH2, totaling ca. 2200 bp. GH genes are easily distinguished by electrophoresis of PCR product, and intron-specific primers exist that preclude cloning (Phillips et al. 2004). These two classes of genes will result in 5kb per specimen for phylogeny reconstruction and AMOVA analyses. Our examination of Mexican trout for these genes indicate a wide range of divergence and are informative in our questions of diversity and evolution. In addition to the Mexican and invasive RT, we will also include in our analysis various divergent taxa within the Cutthroat and Rainbow trout lineages, O. apache, and O. gilae (samples in hand and provided by USFWS). These are the appropriate outgroup taxa following two morphological revisions (Stearley 1992, Stearley & Smith 1993).

Return to Table of Contents

Observed Microsatellite Variation: Preliminary Data—No significant linkage disequilibrium was found among any of the 12 loci examined. All but two populations were in Hardy-Weinberg equilibrium; the two populations deviating from Hardy-Weinberg equilibrium were those from the Granja Truticola hatchery and the wild population immediately below the hatchery on Arroyo la Sidra. Factorial correspondence analysis recovered three groups among the trout populations examined, one from the northern drainages, ríos Yaqui and Casas Grandes, one containing the Mexican golden trout, and one from the southern drainages (ríos Presidio and San Lorenzo) and all hatchery individuals. Variation explained by these groups was low (16.74%), but this was possibly due to smaller sample sizes in the southern drainages adding “noise” to the analysis. AMOVA conducted for these three groups revealed a higher variance between groups than between populations in each group (29% versus 19%). Within population variance comprised 52% of the total variance in the dataset. The high F_st value of 0.4722 indicates significant population substructure. Additional AMOVAs conducted for each of the groups yielded two groups within the Mexican Golden Trout (O. chrysogaster), one composed of trouts from the Río Fuerte and the other contained trout from the ríos Culiacan and Sinaloa. Ten individuals of Mexican Golden Trout from the Río Fuerte exhibited no allelic variation across nine of the twelve loci surveyed. This suggests a recent bottleneck event within this population, with concomitant reductions in genetic variability and potentially lowered levels of fitness. This could also be driving the genetic differentiation within the Río Fuerte and other populations of this species. All pairwise population comparisons using microsatellite markers indicated that hatchery trout were significantly distinct from trout in the adjacent streams. This is an indication that introgression has not completely erased the genetic signature of native trout species. The high F_st identified in the Presidio sample may suggest co-occurrence of native and hatchery fish without introgression or hybridization. Additional specimens from this area are needed to address this question further.

More detailed examination of 545 semi-rare wild trout and cultured Rainbow Trout for two microsatellite loci (Ssa289, One11) again indicated that neither of these loci exhibited linkage disequilibrium, and populations did not deviate from Hardy-Weinberg equilibrium, except hatchery samples and natives found near hatcheries. As with the more complete survey of loci, three genetically homogeneous and distinct groups were found, a northern group (ríos Yaqui-Mayo-Guzmán), Mexican Golden Trout group (Fuerte-Sinaloa-Culiacán) and southern group (ríos Piaxtla-Acaponeta-Baluarte-San Lorenzo). These two loci provided evidence of different degrees of genetic introgression in the genome of the native populations with introduced Rainbow Trout and 2) the persistence of pure stocks of native forms, confirming that native trout are threatened by introgression. We documented hatchery trout in several additional drainages and that introgression is occurring in several locations. The more detailed analysis of the loci Ssa289, One11 almost all major drainages were significantly different from each other; the only exceptions being the ríos Acaponeta-Presidio and the Río San Lorenzo Hatchery–Acaponeta. While this provides weak support for a transfer of fishes between the San Lorenzo hatchery and the Acaponeta, and to a lesser degree into the Presidio, additional loci are required for verification. Notably, F_st and genetic distance between native and hatchery samples from the San Lorenzo (F_st=0.569, D_s=0.499) suggest that native haplotypes in the San Lorenzo have not been completely lost to introgression with hatchery individuals. No evidence of O. mykiss alleles were observed in the Baluarte or Piaxtla drainages; both being rather divergent from all other populations examined, suggesting a long history of isolation.

Genotypic differentiation tests for all pairs of populations revealed that most were significantly different. Those that were not included two populations within the Río Bavispe drainage (northern Río Yaqui), two localities on the Río la Cueva (Río Yaqui) and two different hatchery populations. All pairs of hatchery trout and wild populations immediately adjacent to the hatcheries were significantly different. For these pairs, F_st values were highest for Río Presidio/El Salto hatchery pair (0.45), and lowest for the Río Fuerte/El Aparique hatchery pair (0.05). The Arroyo de la Sidra/Granja Truticola hatchery pair was intermediate with F_st=0.06. One allele of the locus Ssa289 was recovered in Río Yaqui trout but was out of the size class of any of the hatchery fish, indicating a wild origin. Based on the microsatellite data, though genetic introgression of hatchery trout has not obscured the genetic signature of stream trout, the possibility of hybridization cannot be discounted.

While genetic data to date suggest introgression near the ríos Fuerte/El Aparique and Arroyo la Sidra/Granja Truticola hatcheries, they are largely uninformative as to the question of more widespread introgression in these systems. We have collected fish that were morphologically “good” Rainbow Trout in streams in close proximity to both hatcheries, and in sympatry with native trouts. We have also observed early spawning among native populations of trouts.  Specimens collected from the rios Acaponeta, Baluarte, Yaqui and Casas Grandes in early March were in spawning condition; males released milt freely upon removal from the water and females expelled eggs on light pressure. Size structure of these populations indicates that reproductive maturity is reached earlier and at a smaller size than is typical of most populations of native fluvial Oncorhynchus mykiss in the United States. These factors may reduce the impact of hybridization or introgression by reducing temporal overlap in reproductive seasonality between Mexican Oncorhynchus spp. and exotic O. mykiss. This may also produce a skewed directionality of hybridization, resulting from size-assortative pairing, as has been observed in both Rainbow Trout and other species of Oncorhynchus (Foote and Larkin 1988, Dowling and Childs 1992, Rosenfield et al. 2000, Ostberg et al. 2004). While our preliminary data suggest both local hybridization and reductions of variation from bottlenecks, with our sample sizes, additional data are clearly required to accurately identify the long-term impact of O. mykiss-based hatcheries in México.

Return to Table of Contents

Observed DNA Sequence Variation: Preliminary Data—While our analysis did yield limited variation in control region sequences from across much of the natural range of Oncorhynchus mykiss, this marker is very useful for phylogenetic inferences and evaluation of species boundaries.  In our preliminary molecular analysis involving only the control region, 91 out of 1008 characters were parsimony-informative. Mexican Golden Trout (O. chrysogaster) were always embedded within a larger, well-supported clade of samples of the undescribed Mexican trouts. Analysis of partial sequences from O. apache and O. gilae also resulted in those taxa being embedded within these undescribed taxa. Relationships among all forms within the clade inclusive of the Rainbow Trout (O. mykiss) were not resolved in either MP or Bayesian analyses, however, some redband, steelhead and rainbow samples from Bagley and Gall (1999) were consistently resolved with undescribed taxa from the Sierra Madre Occidental in all analyses. Some differentiation among Mexican populations is evident, especially among populations in the ríos Acaponeta/Baluarte/Piaxtla and San Lorenzo/Presidio drainages. This suggests that these populations are native, and have been isolated from northern populations of native (not hatchery) O. mykiss for a considerable time. The ATPase 6/8 sequences also revealed relatively low variation; with only 59 out of 858 parsimony informative characters. Gila Trout and Apache Trout were recovered as sister to the clade containing native (not hatchery) Rainbow Trout. We did recover similar groupings of the southern drainages as found in the control region phylogeny. The nuclear genes GH1C and GH2C were not informative for variation within clade including Rainbow Trout. However, relationships between Gila, Apache, Cutthroat and Rainbow trouts were different from those recovered with mitochondrial genes, especially in the tree recovered using GH1C sequence variation. The genetic variation of these genes across the native trouts of the clade inclusive of native Rainbow Trout and Gila, Apache, and Cutthroat trouts is very promising for ultimately resolving relationships and species boundaries within these clades. Additional genes (ND2, COI, COII, Cyt-b) and microsatellite markers for these and other samples will provide further resolution of the genetic diversity in these river systems.

Figure. Preliminary phylogenetic relationships among populations of undescribed Mexican trout of the "Rainbow Trout Lineage" and the Mexican Golden Trout, O. chrysogaster. Major drainage basins for which samples are available and have been examined are highlighted in different colors. Note that the "Rainbow Trout Lineage" is monophyletic relative to the "Cutthroat Trout Lineage" but that the Mexican Golden Trout is embedded within the taxon currently referred to as O. mykiss ("Rainbow Trout Lineage").

Return to Table of Contents

Phylogeny & Zoogeography 

We will employ a combination of general methods of maximum parsimony and maximum likelihood (Swofford 2000) and Bayesian estimation (implemented with MrBayes v3; Huelsenback & Ronquist 2001) to infer relationships of the taxa. Models of nucleotide substitution will be objectively determined based on likelihood-ratio tests of goodness of fit of all alternative models of sequence evolution (see Sullivan et al. 1997) using an approximation of the null distribution (Yang et al. 1995). Bayesian analyses will be conducted with random starting trees, run 1.0 X 10⁶ generations, and sampled every 100 generations. Base frequencies will be estimated directly from the data whereas initial starting point values for specific nucleotide substitution model parameters will be estimated from the data using an equal-weight parsimony tree (arbitrarily chosen) and specified a priori. In all searches stationarity of the Markov Chain will be determined as the point when sampled log likelihood values plotted against generation time reach a stable mean equilibrium value; "burn-in" data sampled from generations preceding this point were discarded (Huelsenback and Ronquist 2001). All data collected at stationarity will be used to estimate posterior nodal probabilities, mean log likelihood scores, and a summary phylogeny including Bayesian estimates of branch lengths. Five independent replicates will be conducted per data set and inspected for consistency to check for local optima (Huelsenback and Bollback 2001). To ease computational demands, if necessary, we can condense the data matrix to contain only non-redundant haplotypes using COLLAPSE (D. Posada) and perform Bayesian estimations. Nested Clade Analysis will also be conducted for all samples of native trout genotypes. The statistical parsimony algorithm (Posada & Crandall, 2001, Templeton et al. 1992) will be used to evaluate the limits of parsimony with the phylogenetic software TCS (Clement et al. 2000). This can provide an estimate of 95% plausibility set networks within the data set which can then be hierarchically nested and rooted (when appropriate) following the standard methodology (Templeton et al. 1987, Castelloe & Templeton 1994, Crandall 1996, Crandall & Templeton 1993, Templeton & Sing 1993). Interconnections among plausibility networks will then be estimated from a pairwise matrix of absolute number of differences generated with PAUP*v4.0b8 (Swofford 2000). Once all networks are interconnected a separate Bayesian analysis will be conducted (as above) to estimate posterior probabilities of nodes associated with these connections. One can then proceed with the hierarchical nesting algorithm by grouping those clades united with posterior nodal probabilities of 0.95 or higher; nestings poorly supported by both Bayesian posterior probabilities and the statistical parsimony are left unresolved. Significant associations of categorical variables (both sampling localities and morphological diagnostic traits) can then be assessed hierarchically across the nested cladogram using nonparametric, exact permutation contingency tests implemented with GeoDis 2.0 (Posada et al. 2000). For each test, X² will be estimated from a contingency table of categorical data (rows = genetic clades/haplotypes, columns = geographic locations/morphotypes). Significance of observed X² will be determined with an exact test against a simulated null distribution (no association) generated by 10,000 random permutations (Edgington 1987; Roff & Bentzen 1989).

            A more detailed phylogeographic analysis incorporating geographic distances by river miles (km) can be implemented with GeoDis 2.0 (Posada et al. 2000). This interesting approach that can be incorporated with information regarding hybridization and involves the calculation of two distance statistics; clade distance (Dc), measuring geographical spread of a clade (river km), and the nested-clade distance (Dn ), measuring geographic distribution of a clade relative to other clades in it's nesting cohort (Templeton 1995). Interior-tip relationships will be contrasted at each nesting level to approximate the distribution of young versus old haplotypes (Castelloe & Templeton 1994; Crandall & Templeton 1993). Two alternative algorithms for estimating distance statistics are executable with GeoDis 2.0 (Posada et al. 2000); all distance statistics used in this study will be generated using the algorithm which calculates Dc and Dn as unweighted averages from an input file with a pairwise matrix of user-defined population distances. This approach is conceptually consistent with the unweighted methodology outlined by Templeton et al. (1995), whereas the other available algorithm calculates statistics using a weighting scheme based on relative population haploptype/clade frequency. Although the relative performance of these two algorithms has not been thoroughly explored, use of a frequency-based weighting system appears to be sensitive to the effects incomplete and non-uniform sampling (Good & Sullivan 2001). To account for curvature of the earth we will calculate great circle (Templeton 1998), rather than linear, pairwise distances among all populations using the standard great circle formula and Microsoft Excel.

Return to Table of Contents

 Comparative Phylogeography

Quantitative assessment of phylogeographic disjunctions, following the methodologies of Rissler et al 2006), is an effective means of identifying factors driving the diversification of organisms as well as providing a powerful test of alternative biogeographic hypotheses.

Georeferenced data from Miller et al. (2005), our collections to date, and those conducted under the proposed fieldwork will be used to compare the distributions of each species in the study area against its ecological niche model using GARP (Stockwell and Noble 1992, Stockwell and Peters 1999). Individuals sampled for phylogenetic analysis (for trout, Campostoma, Codoma, Gila spp., and Catostomus spp.) will be used to identify the geographic distribution of lineages within each taxon, and (when lineages have allopatric distributions) used to infer lineage membership among the larger data set of observations. As these are aquatic organisms and more strongly tied to hydrological features, we expect that the spatial distribution of breaks will be more fine-grained than those observed in terrestrial amphibians by Rissler et al (2006); we will thus use a finer grid than those authors (10x10 km vs 25x25km) to identify geographic breaks between lineages. Areas with a higher than random level of phylogeographic breaks will be identified by comparing the number of breaks in each cell against an expected value (number of breaks in study area/number of species in study area). These areas will then be scrutinized for possible geological or geographic barriers. Given the high rate of diversification by vicariance in freshwater stream fishes, we suspect a strong pattern of congruence among breaks.

Return to Table of Contents

Amphibian Biodiversity in Sierra Madre Occidental

Return to Table of Contents