BACKGROUND
originating independently in diverse taxa, the xx/xy sex determination system is a striking example of convergent evolution  <cit> .
in many species with y chromosomes, a y-linked dominant male-determining gene is hypothesized to initiate male sexual differentiation  <cit> .
however, y-linked dominant male-determining genes have eluded characterization with only a few exceptions.
in therian mammals the dominant male-determining factor, sry, is one of the few y-linked genes retained after the degeneration of the y chromosome  <cit> .
in the medaka fish oryzias latipes the y-linked gene dmy functions as the dominant male determining gene  <cit> .
numerous insect species determine sex through the xx/xy sex-determination system, although insect y chromosomes likely originated independently multiple times from diverse origins  <cit> .
a notable exception is that while drosophila species have y chromosomes, sex is not determined by a gene on the y, but instead the x-to-autosome ratio  <cit> .
although no y-linked sex-determining genes have been characterized in insects, y-linked sex-determining genes likely exist in anopheles mosquitoes, the housefly musca domestica, the mediterranean fruit fly ceratitis capitata, and possibly many more species with y chromosomes  <cit> .
most insect y chromosomes remain entirely unexplored.
over time, male-beneficial genes may accumulate on the y chromosome by gene duplications  <cit> .
for example, all known protein-coding genes on the d. melanogaster y chromosome originated from autosome-to-y duplications  <cit> .
due to their exclusive presence in males, y chromosomes may also serve as a sanctuary for genes advantageous to males, and deleterious to females  <cit> .
y chromosome genes have eluded all but the most determined attempts at characterization due to the heterochromatic and repetitive nature of the y  <cit> .
using bacterial artificial chromosomes  and an iterative mapping and sequencing strategy, the sequence of the euchromatic region of the y chromosome was obtained for humans and two other primates  <cit> .
the iterative mapping strategy used in primates has not been implemented to sequence other y chromosomes due to the cost associated with the technique.
in insects, the sequences of y genes have only been characterized in a few drosophila species.
the drosophila y gene sequences were found by a strategy, developed by carvalho and colleagues, based on the comparison of sequences not anchored to chromosomes to known protein sequences  <cit> .
this method may be difficult to implement in species without a reference genome where most of the sequences are anchored to chromosomes, and it failed to identify any y genes in the african malaria mosquito an.
gambiae <cit> .
in this study, we introduce the chromosome quotient  method, a novel approach to systematically discover y genes.
the cq method takes advantage of high-throughput genome and transcriptome sequencing data and does not require a reference genome where most of the sequences are anchored to chromosomes.
we used the cq method to find y sequences from two important malaria mosquitoes, an.
stephensi and an.
gambiae, members of the same subgenus cellia <cit> .
comparison of the y sequences to transcriptome sequence data led to the discovery and characterization of six anopheles y genes, three from each species.
evolutionary analysis strongly suggests rapid y chromosome evolution in anopheles mosquitoes.
RESULTS
the chromosome quotient  method
we developed the cq method based on the following observations.
y chromosomes are present only in males.
therefore, sequences unique to the y should be present in male sequence data and absent from female sequence data.
using this principal, we developed a method to identify y sequences by comparing separate male and female sequence data to a reference genome.
unique y sequences can be identified because they are only present in the male sequence data.
however, searching for sequences exclusive to the male sequence data yielded only a few short sequences because most y sequences appear to be at least partially repetitive.
we realized that if the reference genome was fragmented into smaller pieces by removing repeats, more sequences exclusive to the male sequence data could be identified.
y chromosomes appear to be primarily composed of repetitive sequences that are not exclusive to the y. since most y-repeats have closely-related sequences on the autosome or x, they may appear be present in both the female and male sequence data.
recent duplications of autosome or x sequences to the y may also cause y sequences to appear to be present in both the male and female sequence data.
we reduced the interference from repetitive sequences by using very strict alignment criteria.
we require zero mismatches over the entire length of the reads from the male and female sequence data.
we increased the number of y chromosome sequences identified by allowing a few alignments from the female sequence data, as long as there were many more alignments from the male sequence data.
to define the number of alignments allowed from female data we use the ratio of female to male alignments, a parameter we call the chromosome quotient .
for a given sequence si, cq = f/m, where f is the number of alignments from female sequence data to si, and m is the number of alignments from male sequence data to si.
the method by which chromosome quotients are calculated is found in the additional files .
the chromosome quotient allows for the differentiation of y sequences from autosome and x sequences .
males and females share the same complement of autosomes, so autosomal sequences are present in both male and female sequence data in roughly the same quantities.
therefore, autosomal sequences have chromosome quotients distributed around one.
females have two x chromosomes while males have only one, so x chromosome sequences are present in female sequence data roughly twice as frequently as in male sequence data.
therefore, x chromosome sequences have chromosome quotients distributed around two.
unique y sequences are present only in male sequence data, and therefore have chromosome quotients of zero.
repetitive y sequences are present in both the male and female sequence data and have chromosome quotients greater than zero but less than the chromosome quotients of autosomal sequences.
we set a threshold of cq =  <dig>  to differentiate y sequences from autosome and x sequences.
using the cq method, a sequence is classified as y-linked if it has more than  <dig>  times as many alignments from male data than from female data, and therefore a cq less than  <dig> .
in case the coverage of the male and female sequence data differs, chromosome quotients are normalized to the median chromosome quotient of known autosomal sequences.
y chromosome sequences have distinctive chromosome quotients
to validate that the chromosome quotients of y sequences are distinctive from those of autosome and x sequences, the cq method was tested in h. sapiens and d. melanogaster.
autosome, x, and y sequences were downloaded for both h. sapiens and d. melanogaster.
the repetitive sequences indicated by repeatmasker were removed, fragmenting the genomes into many smaller pieces.
male and female sequence data was located for both species.
the chromosome quotients of the fragmented autosome, x, and y sequences were calculated using the separate male and female sequence data .
the number of sequences and the normalized median chromosome quotients  of the positive control data sets illustrate that autosomal sequences have chromosome quotients distributed around one, x sequences have chromosome quotients distributed around two, and y sequences have chromosome quotients near zero.
the distributions of chromosome quotients from different chromosomes are shown in figure  <dig>
in both species, autosomal sequences have chromosome quotients distributed around one, x sequences have chromosome quotients distributed around two, and most y sequences have chromosome quotients near zero.
the median chromosome quotients of autosome and x sequences are significantly greater than the median chromosome quotients of y sequences .
thus, y sequences can be identified by their distinctive chromosome quotients.
the interval of chromosome quotients  captures more than  <dig> percent of the known y sequences with a rate of false positives less than  <dig>  percent in h. sapiens and d. melanogaster .
therefore, we chose this threshold for the discovery of y chromosome sequences in the mosquitoes an.
stephensi and an.
gambiae.
the false positive rate is highly dependent on the cq cut-off below which a sequence is classified as y-linked.
in this study, we used a cq cut-off of  <dig>  as it minimized false positives while keeping false negatives relatively low.
there is flexibility in selecting the cq cut-off according to the focus of the study and the scale of downstream verification that the investigators plan to perform.
discovery of novel y sequences in anopheles mosquitoes
we used the cq method to discover novel y sequences in two malaria mosquitoes.
in this study, we consider a sequence to be y-linked if it has a chromosome quotient less than  <dig> , meaning it has more than  <dig>  times more alignments from male data than from female data.
we first searched for y sequences in the asian malaria mosquito an.
stephensi, a species with a draft genome generated by  <dig> sequencing .
we performed illumina sequencing on male and female an.
stephensi genomic dna .
chromosome quotients were first calculated for all known autosome and x contigs in the an.
stephensi genome .
all of the autosome and x sequences had chromosome quotients greater than  <dig> .
chromosome quotients were then calculated for  <dig>  contigs from the an.
stephensi genome.
we identified  <dig> candidate y sequences with chromosome quotients less than  <dig> , encompassing more than  <dig>  bases .
five of these sequences were tested for male-specific amplification with pcr on male and virgin-female genomic dna .
all five sequences amplified a product exclusively in male genomic dna .
all of the sequences listed amplified a pcr product exclusively in male genomic dna, and all have chromosome quotients less than  <dig> .
we made an exception for syg <dig> as it only has  <dig> alignments from male data, and the cutoff was  <dig> .
even at lower stringency, syg <dig> still had no female alignments so we decided to include it.
the chromosome quotients for gyg <dig> and gyg <dig> were calculated from the transcripts because of the long introns present in the genes.
we then applied the cq method to the genome of the african malaria mosquito, an.
gambiae.
we performed illumina sequencing on male and female an.
gambiae genomic dna .
the an.
gambiae agamp <dig> genome assembly was downloaded from vectorbase.
repetitive sequences identified by repeatmasker were removed, fragmenting the genome into many smaller pieces.
chromosome quotients were calculated for the fragmented autosome and x sequences from the agamp <dig> genome assembly .
the autosomal sequences had chromosome quotients distributed around one, and the x chromosome sequences had chromosome quotients distributed around two.
chromosome quotients were calculated for the fragmented y chromosome sequences from the agamp <dig> genome assembly.
of the  <dig> known an.
gambiae y sequences,  <dig> had chromosome quotients less than  <dig> .
the an.
gambiae sequences that were unable to be anchored to chromosomes, known as unkn in the agamp <dig> assembly, were fragmented by removing the repeats indicated by repeatmasker.
chromosome quotients were calculated for these unkn sequences.
from the unkn sequences, we identified  <dig> novel candidate y sequences, encompassing more than  <dig>  bases .
we also tested the cq method de novo by generating a genome assembly from male g <dig> strain illumina sequence data using abyss.
chromosome quotients were calculated for all the sequences in the illumina assembly.
from the illumina assembly, we identified  <dig> novel candidate y sequences encompassing more than  <dig>  bases .
all y sequences reported in this study are provided in the additional files .
discovery and expression profile of six anopheles y genes
comparison of the y sequences to an.
stephensi and an.
gambiae rna-sequencing data with blastn led to the identification of six novel y genes, three from each species .
the an.
stephensi and an.
gambiae y genes were given the names syg  and gyg  respectively, followed by a number designating the order in which they were discovered.
all six y genes amplified a product exclusively in male genomic dna .
sequencing of rt-pcr products  confirmed transcription of all six genes.
cdna sequences for genes with confirmed introns were used to determine the structure of these genes, and are provided in additional file  <dig> .
interestingly, all six genes were expressed in the early embryo, 2–12 hours after egg deposition .
none of the six genes were expressed in 0–1 hour embryos, when only maternal transcripts are present, or in the adult female where no y chromosome is present  <cit> .
the three an.
stephensi y genes include two recent duplications to the y and a novel y gene
syg <dig>  shares  <dig> percent nucleotide identity with a gene mapped to chromosome 2r .
we call this autosomal paralog of syg <dig>  asa-bbx .
both syg <dig> and asa-bbx are homologous to the d. melanogaster gene bobby sox  and to two closely linked an.
gambiae genes, agap <dig> and agap <dig>  which are located on the an.
gambiae 2r chromosome arm.
asa-bbx is part of a large contig and the synteny of its neighboring genes is conserved among anopheles and aedes mosquitoes  suggesting that asa-bbx is the ancestral copy, and syg <dig> was formed from a recent duplication of asa-bbx to the y. phylogenetic analysis suggests that this duplication of syg <dig> occurred after the evolutionary divergence of an.
stephensi and an.
gambiae .
5′ and 3′ race was performed on syg <dig> to characterize its full-length cdna .
copy number determination using digital droplet pcr  <cit> , designed to detect both syg <dig> and asa-bbx, indicates two autosomal asa-bbx copies per diploid genome  and six-to-nine copies of syg <dig> in males .
chromosomal fluorescence in situ hybridization  with probes designed according to syg <dig> hybridized to the an.
stephensi y on male mitotic chromosomes and on the 2r arm on female polytene chromosomes .
similar to the digital pcr analysis, this fish result confirms that asa-bbx is a single-copy autosomal gene located in subdivision 18a of the 2r arm, while syg <dig> is a multi-copy gene on the y.
the most similar non-y sequences to the y genes in the an.
stephensi and an.
gambiae genomes.
<dig>  query sequences are the nucleotide sequences of the y genes shown in additional file  <dig>   <dig>  blastn is used to obtain e-values.
percent identities reflect the entire sequence except for gyg <dig> and gyg <dig>   <dig>  asa-bbx is orthologous to an.
gambiae agap003896/ <dig>  asa-bbx is an autosomal paralog of syg <dig>   <dig>  syg <dig> has no similar sequence in the rest of the an.
stephensi genome except for a 55 bp repetitive region.
<dig>  syg <dig> is closely related to an annotated gene in an.
stephensi, which is orthologous to an.
gambiae agap <dig>   <dig>  gyg <dig> appears to be a composite of different fragments that share high similarities to different un-annotated non-y sequences in an.
gambiae.
some fragments are repeated in the genome.
shown in this table is a  <dig> base fragment that is  <dig> percent identical to sequences in chromosome 2r, 3 l as well as unmapped scaffolds.
<dig>  gyg <dig> and gyg <dig> overlap, which is why the best non-y match is identical between the two y genes.
gyg <dig> and gyg <dig> are transcribed from opposite strands.
<dig>  gyg <dig> is most similar to an un-annotated 2 l sequence and an unmapped sequence in an.
gambiae.
with the exception of a  <dig> base repetitive region, syg <dig>  matched no other sequences in the an.
stephensi genome using blastn or tblastx with an e-value threshold of 1e- <dig>  syg <dig> has no homolog in an.
gambiae as indicated by blastn and tblastx searches of the pest, m, and s genome assemblies as well as the trace files.
syg <dig>  shares  <dig>  percent nucleotide identity with an autosomal or x gene within the  <dig>  base scaffold <dig> of the an.
stephensi genome.
this gene is homologous to agap <dig> on the an.
gambiae x chromosome .
agap <dig> and its neighboring genes showed conserved synteny in an.
gambiae, an.
stephensi, and culex quinquefasciatus, suggesting that the autosomal or x paralog of syg <dig> is ancestral.
thus, we deduce syg <dig> is also a recent duplication to the an.
stephensi y that occurred after the divergence of an.
stephensi and an.
gambiae.
the annotated homologs of the y genes provide clues to putative functions of the y genes we discovered.
<dig>  homologs were identified using blastx against the an.
gambiae protein database  and the ncbi non-redundant protein database.
not surprisingly, the best matches are to an.
gambiae sequences because the queries are from an.
gambiae and an.
stephensi.
the an.
gambiae genes agap <dig> and agap <dig> are orthologous to the autosomal or x paralogs of syg <dig> and syg <dig>  respectively.
the other homologs listed in column  <dig> are either distantly related or only share partial overlap.
in cases where the y genes match more than one related homolog, the homolog with the best e-value is shown.
the three homologs shown for gyg <dig> are unrelated.
<dig>  maximum percent identities and e-values were obtained by blastx against the an.
gambiae protein database .
for y genes that have confirmed introns , transcript sequences were used as the query.
<dig>  the match to agap <dig> is very short at the amino acid level although the nucleotide identity is  <dig> percent, and  <dig> percent in two longer fragments .
the sense strand of gyg <dig> is reverse-complementary to agap <dig>   <dig>  the sense strand of gyg <dig> is reverse-complementary to agap <dig>
the three an.
gambiae y genes have closely related autosomal or x paralogs in the an.
gambiae genome
gyg <dig>  and gyg <dig>  are two new genes found in the y chromosome sequences of the agamp <dig> genome assembly.
although they overlap on opposite strands, there are unique transcribed regions that allowed specific rt-pcr products to be amplified for each gene .
both genes are highly transcribed from the early embryo until the adult male, including in the testes.
sequencing of race and rt-pcr products confirmed the organization of these two genes .
both genes have closely-related autosomal or x paralogs although these paralogs are not annotated .
homologs of gyg <dig> and gyg <dig> were not identified in an.
stephensi by blastn and tblastx searches of the genome assemblies as well as the  <dig> trace files, although a short repetitive region of gyg <dig> does have significant alignments to the an.
stephensi genome from tblastx.
a portion of the gyg <dig> transcript shares more than  <dig> percent nucleotide identity with the agap <dig> group of genes on both the x and autosomes, while the flanking regions are distantly but significantly related to agap <dig> and agap <dig>  respectively.
a portion of the gyg <dig> transcript shares greater than  <dig> percent nucleotide identity with a multi-copy non-coding sequence in the an.
gambiae genome.
another part of gyg <dig> also contains a region that is homologous to agap <dig> .
the gyg <dig> and gyg <dig> transcripts are reverse complementary to agap <dig> and agap <dig>  respectively.
this transcription orientation coupled with the fact that gyg <dig> and gyg <dig> are in opposite orientation may indicate functions related to double-stranded rna formation, rather than coding for proteins.
because there are no known orthologs to agap <dig>  agap <dig>  and agap <dig> we could not determine the ancestry of the gyg <dig> and gyg <dig> genes relative to their paralogs.
gyg <dig>  was discovered in an assembly of our own an.
gambiae illumina data.
there are two similar sequences in the an.
gambiae genome, both sharing  <dig> percent nucleotide identity with gyg <dig> .
one is on chromosome 2 l and the other is not mapped but its chromosome quotient suggests it is an autosomal or x sequence.
thus, gyg <dig> is again involved in a recent duplication, although we cannot assign the direction of the duplication.
a blastx search shows that gyg <dig> is distantly related to agap <dig> , which belongs to a large family of general transcription factor ii-i repeat domain-containing proteins .
discussion
here we have shown that the cq method can reliably identify y sequences using the ratio of alignments from male and female sequence data.
using the cq method we identified  <dig> y chromosome sequences encompassing more than  <dig>  bases of novel y chromosome sequences from the malaria mosquitoes an.
stephensi and an.
gambiae.
it is also possible to identify y sequences by coverage depth  <cit> .
in male sequence data, y sequences should have half the coverage of autosomal sequences.
in mixed male and female sequence data, y sequences should have one-quarter the coverage of autosomal sequences  <cit> .
the cq method does not directly rely on coverage depth to identify y sequences.
while y sequences often do have low coverage depth, it is often not sufficient to classify a sequence as y-linked solely by low coverage depth.
mis-assemblies of the reference genome, and allelic differences between the reference genome and the short-read sequence data could lead to low coverage depth causing errant classification of y-linkage.
unmasked repetitive sequences may lead to high coverage depth, obfuscating many y sequences.
instead, the cq method uses the ratio of female to male alignments to a reference sequence.
the cq method takes into account both the lack of female alignments and the presence of male alignments when classifying sequences as y-linked, reducing the overall rate of false positives and false negatives.
mis-assemblies and allelic variation can be detected by the lack of male reads aligned to a sequence.
for a sequence to be classified as y-linked with the chromosome quotient it must have many alignments from males, and zero or a few alignments from females.
although the functions of the y genes that we found are not yet clear, there are some tantalizing hints to their potential functions.
syg <dig> is homologous to an hmg-box transcription factor involved in a number of developmental functions.
hmg-box transcription factors include the two key genes involved in male-determination in mammals  <cit> .
syg <dig> shares homology with a general transcription factor ii repeat domain containing protein, indicating the possibility that it is a y-linked transcription factor.
five of the y genes are clearly expressed in the testes suggesting that they may perform male-specific roles.
all six y genes are expressed in the early embryo, which is the time when sexual differentiation begins  <cit> .
our research offers an opportunity to investigate y chromosome evolution, which until now has been restricted to limited lineages.
primate y chromosomes have been shown to undergo strict conservation following rapid gene loss  <cit> .
in d. melanogaster, however, almost all the y genes are the result of duplications to the y. here we have shown that two of the three an.
stephensi y genes were recently acquired from the autosomes or the x. all three an.
gambiae y genes have highly similar paralogs on the autosomes or x , again consistent with recent duplication.
thus, gene acquisitions through duplications appear to be a major source of y genes in anopheles mosquitoes.
this is intriguing given that unlike d. melanogaster, the anopheles y is thought to harbor a dominant male-determining factor that initiates sexual differentiation, similar to humans and several non-drosophila flies  <cit> .
so far, no conservation between the genes on the two anopheles y chromosomes has been found, although the two species belong to the same subgenus.
it is possible that we simply missed the common y genes that originated prior to the divergence between an.
gambiae and an.
stephensi.
however, such shared y genes are less likely to have highly similar autosomal or x paralogs compared to recent duplications, and they may even be unique to the y. y genes that are unique or significantly diverged from their autosomal or x paralogs are the easiest to find with the cq method.
to mitigate the effect of potential underrepresentation of y genes in the anopheles genome sequences, we compared the trace sequence files of both species and found no indications of shared y genes.
thus, the evidence so far suggests that the an.
gambiae and an.
stephensi y chromosomes are very different in their gene content.
this conclusion is consistent with the observation that two of the three y genes in an.
stephensi were recently acquired after the divergence of the two mosquitoes.
further investigation of y genes in additional anopheles species will reveal whether the evolutionary scenarios demonstrated by the drosophila y chromosomes , namely recent y gene acquisition and y chromosome replacement, are widely applicable in species with a dominant male-determining factor on the y.
some of the y genes described in this study have unique y chromosome sequences.
because of their exclusive paternal inheritance and general lack of recombination, these y sequences could serve as powerful molecular markers for the investigation of population structure and incipient speciation, which are often relevant to vectorial capacity  <cit> .
only female mosquitoes bite and transmit pathogens.
genetic strategies to control mosquito-borne diseases may include the release of sterile males , the release of insects carrying a dominant lethal gene , or the replacement of a pathogen-susceptible population with a pathogen-resistant one  <cit> .
under these scenarios, the release of males is either required or preferred.
a better understanding of the possible functions of mosquito y genes in sexual differentiation will facilitate the production of sterile males for improved sit.
it may also enable the production of highly competitive males to improve sit and ridl and to help spread pathogen resistance genes.
CONCLUSIONS
in this study, we have shown that y sequences can be identified by their distinctive chromosome quotients, the ratio of female to male alignments.
although the genomes of many organisms with y chromosomes have been sequenced, y sequences have been characterized from only a small subset of these genomes.
y sequences are surely present in these genomes, but are probably fragmented into short contigs challenging the traditional experimental methods of finding y sequences.
the quality of the genome assembly is not a limiting factor in finding y sequences, making the cq method widely applicable to species with fragmented genome.
autosome or x sequences rarely have chromosome quotients less than  <dig> , giving the cq method a low false positive rate.
the cq method is able to identify the majority of known y sequences from species we have analyzed.
we have also shown that the cq method can identify recent duplications to the y that share more than  <dig> percent nucleotide identity with autosomal or x paralogs.
coupling the cq method with transcriptome sequencing allowed us to identify six y genes in an.
stephensi and an.
gambiae.
analysis of the six y genes indicates rapid y chromosome evolution between an.
stephensi and an.
gambiae.
the anopheles y genes characterized in this study provide unique markers for population and phylogenetic analysis, and opportunities for novel mosquito control measures through the manipulation of sexual dimorphism and fertility.
the identification of y genes in more organisms will help us understand the evolution and functions of y chromosomes.

 originating independently in diverse taxa, the xx/xy sex determination system is a striking example of convergent evolution. in many species with y chromosomes, a y-linked dominant male-determining gene is hypothesized to initiate male sexual differentiation. however, y-linked dominant male-determining genes have eluded characterization with only a few exceptions. in therian mammals the dominant male-determining factor, sry, is one of the few y-linked genes retained after the degeneration of the y chromosome. in the medaka fish oryzias latipes the y-linked gene dmy functions as the dominant male determining gene. numerous insect species determine sex through the xx/xy sex-determination system, although insect y chromosomes likely originated independently multiple times from diverse origins. a notable exception is that while drosophila species have y chromosomes, sex is not determined by a gene on the y, but instead the x-to-autosome ratio. although no y-linked sex-determining genes have been characterized in insects, y-linked sex-determining genes likely exist in anopheles mosquitoes, the housefly musca domestica, the mediterranean fruit fly ceratitis capitata, and possibly many more species with y chromosomes. most insect y chromosomes remain entirely unexplored. over time, male-beneficial genes may accumulate on the y chromosome by gene duplications. for example, all known protein-coding genes on the d. melanogaster y chromosome originated from autosome-to-y duplications. due to their exclusive presence in males, y chromosomes may also serve as a sanctuary for genes advantageous to males, and deleterious to females. y chromosome genes have eluded all but the most determined attempts at characterization due to the heterochromatic and repetitive nature of the y. using bacterial artificial chromosomes and an iterative mapping and sequencing strategy, the sequence of the euchromatic region of the y chromosome was obtained for humans and two other primates. the iterative mapping strategy used in primates has not been implemented to sequence other y chromosomes due to the cost associated with the technique. in insects, the sequences of y genes have only been characterized in a few drosophila species. the drosophila y gene sequences were found by a strategy, developed by carvalho and colleagues, based on the comparison of sequences not anchored to chromosomes to known protein sequences. this method may be difficult to implement in species without a reference genome where most of the sequences are anchored to chromosomes, and it failed to identify any y genes in the african malaria mosquito an. gambiae. in this study, we introduce the chromosome quotient method, a novel approach to systematically discover y genes. the cq method takes advantage of high-throughput genome and transcriptome sequencing data and does not require a reference genome where most of the sequences are anchored to chromosomes. we used the cq method to find y sequences from two important malaria mosquitoes, an. stephensi and an. gambiae, members of the same subgenus cellia. comparison of the y sequences to transcriptome sequence data led to the discovery and characterization of six anopheles y genes, three from each species. evolutionary analysis strongly suggests rapid y chromosome evolution in anopheles mosquitoes.  the chromosome quotient method we developed the cq method based on the following observations. y chromosomes are present only in males. therefore, sequences unique to the y should be present in male sequence data and absent from female sequence data. using this principal, we developed a method to identify y sequences by comparing separate male and female sequence data to a reference genome. unique y sequences can be identified because they are only present in the male sequence data. however, searching for sequences exclusive to the male sequence data yielded only a few short sequences because most y sequences appear to be at least partially repetitive. we realized that if the reference genome was fragmented into smaller pieces by removing repeats, more sequences exclusive to the male sequence data could be identified. y chromosomes appear to be primarily composed of repetitive sequences that are not exclusive to the y. since most y-repeats have closely-related sequences on the autosome or x, they may appear be present in both the female and male sequence data. recent duplications of autosome or x sequences to the y may also cause y sequences to appear to be present in both the male and female sequence data. we reduced the interference from repetitive sequences by using very strict alignment criteria. we require zero mismatches over the entire length of the reads from the male and female sequence data. we increased the number of y chromosome sequences identified by allowing a few alignments from the female sequence data, as long as there were many more alignments from the male sequence data. to define the number of alignments allowed from female data we use the ratio of female to male alignments, a parameter we call the chromosome quotient. for a given sequence si, cq = f/m, where f is the number of alignments from female sequence data to si, and m is the number of alignments from male sequence data to si. the method by which chromosome quotients are calculated is found in the additional files. the chromosome quotient allows for the differentiation of y sequences from autosome and x sequences. males and females share the same complement of autosomes, so autosomal sequences are present in both male and female sequence data in roughly the same quantities. therefore, autosomal sequences have chromosome quotients distributed around one. females have two x chromosomes while males have only one, so x chromosome sequences are present in female sequence data roughly twice as frequently as in male sequence data. therefore, x chromosome sequences have chromosome quotients distributed around two. unique y sequences are present only in male sequence data, and therefore have chromosome quotients of zero. repetitive y sequences are present in both the male and female sequence data and have chromosome quotients greater than zero but less than the chromosome quotients of autosomal sequences. we set a threshold of cq =  to differentiate y sequences from autosome and x sequences. using the cq method, a sequence is classified as y-linked if it has more than times as many alignments from male data than from female data, and therefore a cq less than. in case the coverage of the male and female sequence data differs, chromosome quotients are normalized to the median chromosome quotient of known autosomal sequences. y chromosome sequences have distinctive chromosome quotients to validate that the chromosome quotients of y sequences are distinctive from those of autosome and x sequences, the cq method was tested in h. sapiens and d. melanogaster. autosome, x, and y sequences were downloaded for both h. sapiens and d. melanogaster. the repetitive sequences indicated by repeatmasker were removed, fragmenting the genomes into many smaller pieces. male and female sequence data was located for both species. the chromosome quotients of the fragmented autosome, x, and y sequences were calculated using the separate male and female sequence data. the number of sequences and the normalized median chromosome quotients of the positive control data sets illustrate that autosomal sequences have chromosome quotients distributed around one, x sequences have chromosome quotients distributed around two, and y sequences have chromosome quotients near zero. the distributions of chromosome quotients from different chromosomes are shown in figure  in both species, autosomal sequences have chromosome quotients distributed around one, x sequences have chromosome quotients distributed around two, and most y sequences have chromosome quotients near zero. the median chromosome quotients of autosome and x sequences are significantly greater than the median chromosome quotients of y sequences. thus, y sequences can be identified by their distinctive chromosome quotients. the interval of chromosome quotients captures more than percent of the known y sequences with a rate of false positives less than percent in h. sapiens and d. melanogaster. therefore, we chose this threshold for the discovery of y chromosome sequences in the mosquitoes an. stephensi and an. gambiae. the false positive rate is highly dependent on the cq cut-off below which a sequence is classified as y-linked. in this study, we used a cq cut-off of as it minimized false positives while keeping false negatives relatively low. there is flexibility in selecting the cq cut-off according to the focus of the study and the scale of downstream verification that the investigators plan to perform. discovery of novel y sequences in anopheles mosquitoes we used the cq method to discover novel y sequences in two malaria mosquitoes. in this study, we consider a sequence to be y-linked if it has a chromosome quotient less than, meaning it has more than times more alignments from male data than from female data. we first searched for y sequences in the asian malaria mosquito an. stephensi, a species with a draft genome generated by sequencing. we performed illumina sequencing on male and female an. stephensi genomic dna. chromosome quotients were first calculated for all known autosome and x contigs in the an. stephensi genome. all of the autosome and x sequences had chromosome quotients greater than. chromosome quotients were then calculated for contigs from the an. stephensi genome. we identified candidate y sequences with chromosome quotients less than, encompassing more than bases. five of these sequences were tested for male-specific amplification with pcr on male and virgin-female genomic dna. all five sequences amplified a product exclusively in male genomic dna. all of the sequences listed amplified a pcr product exclusively in male genomic dna, and all have chromosome quotients less than. we made an exception for syg as it only has alignments from male data, and the cutoff was. even at lower stringency, syg still had no female alignments so we decided to include it. the chromosome quotients for gyg and gyg were calculated from the transcripts because of the long introns present in the genes. we then applied the cq method to the genome of the african malaria mosquito, an. gambiae. we performed illumina sequencing on male and female an. gambiae genomic dna. the an. gambiae agamp genome assembly was downloaded from vectorbase. repetitive sequences identified by repeatmasker were removed, fragmenting the genome into many smaller pieces. chromosome quotients were calculated for the fragmented autosome and x sequences from the agamp genome assembly. the autosomal sequences had chromosome quotients distributed around one, and the x chromosome sequences had chromosome quotients distributed around two. chromosome quotients were calculated for the fragmented y chromosome sequences from the agamp genome assembly. of the known an. gambiae y sequences, had chromosome quotients less than. the an. gambiae sequences that were unable to be anchored to chromosomes, known as unkn in the agamp assembly, were fragmented by removing the repeats indicated by repeatmasker. chromosome quotients were calculated for these unkn sequences. from the unkn sequences, we identified novel candidate y sequences, encompassing more than bases. we also tested the cq method de novo by generating a genome assembly from male g strain illumina sequence data using abyss. chromosome quotients were calculated for all the sequences in the illumina assembly. from the illumina assembly, we identified novel candidate y sequences encompassing more than bases. all y sequences reported in this study are provided in the additional files. discovery and expression profile of six anopheles y genes comparison of the y sequences to an. stephensi and an. gambiae rna-sequencing data with blastn led to the identification of six novel y genes, three from each species. the an. stephensi and an. gambiae y genes were given the names syg and gyg respectively, followed by a number designating the order in which they were discovered. all six y genes amplified a product exclusively in male genomic dna. sequencing of rt-pcr products confirmed transcription of all six genes. cdna sequences for genes with confirmed introns were used to determine the structure of these genes, and are provided in additional file. interestingly, all six genes were expressed in the early embryo, 2–12 hours after egg deposition. none of the six genes were expressed in 0–1 hour embryos, when only maternal transcripts are present, or in the adult female where no y chromosome is present. the three an. stephensi y genes include two recent duplications to the y and a novel y gene syg shares percent nucleotide identity with a gene mapped to chromosome 2r. we call this autosomal paralog of syg asa-bbx. both syg and asa-bbx are homologous to the d. melanogaster gene bobby sox and to two closely linked an. gambiae genes, agap and agap which are located on the an. gambiae 2r chromosome arm. asa-bbx is part of a large contig and the synteny of its neighboring genes is conserved among anopheles and aedes mosquitoes suggesting that asa-bbx is the ancestral copy, and syg was formed from a recent duplication of asa-bbx to the y. phylogenetic analysis suggests that this duplication of syg occurred after the evolutionary divergence of an. stephensi and an. gambiae. 5′ and 3′ race was performed on syg to characterize its full-length cdna. copy number determination using digital droplet pcr, designed to detect both syg and asa-bbx, indicates two autosomal asa-bbx copies per diploid genome and six-to-nine copies of syg in males. chromosomal fluorescence in situ hybridization with probes designed according to syg hybridized to the an. stephensi y on male mitotic chromosomes and on the 2r arm on female polytene chromosomes. similar to the digital pcr analysis, this fish result confirms that asa-bbx is a single-copy autosomal gene located in subdivision 18a of the 2r arm, while syg is a multi-copy gene on the y. the most similar non-y sequences to the y genes in the an. stephensi and an. gambiae genomes. query sequences are the nucleotide sequences of the y genes shown in additional file blastn is used to obtain e-values. percent identities reflect the entire sequence except for gyg and gyg asa-bbx is orthologous to an. gambiae agap003896/ asa-bbx is an autosomal paralog of syg syg has no similar sequence in the rest of the an. stephensi genome except for a 55 bp repetitive region. syg is closely related to an annotated gene in an. stephensi, which is orthologous to an. gambiae agap gyg appears to be a composite of different fragments that share high similarities to different un-annotated non-y sequences in an. gambiae. some fragments are repeated in the genome. shown in this table is a base fragment that is percent identical to sequences in chromosome 2r, 3 l as well as unmapped scaffolds. gyg and gyg overlap, which is why the best non-y match is identical between the two y genes. gyg and gyg are transcribed from opposite strands. gyg is most similar to an un-annotated 2 l sequence and an unmapped sequence in an. gambiae. with the exception of a base repetitive region, syg matched no other sequences in the an. stephensi genome using blastn or tblastx with an e-value threshold of 1e- syg has no homolog in an. gambiae as indicated by blastn and tblastx searches of the pest, m, and s genome assemblies as well as the trace files. syg shares percent nucleotide identity with an autosomal or x gene within the base scaffold of the an. stephensi genome. this gene is homologous to agap on the an. gambiae x chromosome. agap and its neighboring genes showed conserved synteny in an. gambiae, an. stephensi, and culex quinquefasciatus, suggesting that the autosomal or x paralog of syg is ancestral. thus, we deduce syg is also a recent duplication to the an. stephensi y that occurred after the divergence of an. stephensi and an. gambiae. the annotated homologs of the y genes provide clues to putative functions of the y genes we discovered. homologs were identified using blastx against the an. gambiae protein database and the ncbi non-redundant protein database. not surprisingly, the best matches are to an. gambiae sequences because the queries are from an. gambiae and an. stephensi. the an. gambiae genes agap and agap are orthologous to the autosomal or x paralogs of syg and syg respectively. the other homologs listed in column are either distantly related or only share partial overlap. in cases where the y genes match more than one related homolog, the homolog with the best e-value is shown. the three homologs shown for gyg are unrelated. maximum percent identities and e-values were obtained by blastx against the an. gambiae protein database. for y genes that have confirmed introns, transcript sequences were used as the query. the match to agap is very short at the amino acid level although the nucleotide identity is percent, and percent in two longer fragments. the sense strand of gyg is reverse-complementary to agap the sense strand of gyg is reverse-complementary to agap the three an. gambiae y genes have closely related autosomal or x paralogs in the an. gambiae genome gyg and gyg are two new genes found in the y chromosome sequences of the agamp genome assembly. although they overlap on opposite strands, there are unique transcribed regions that allowed specific rt-pcr products to be amplified for each gene. both genes are highly transcribed from the early embryo until the adult male, including in the testes. sequencing of race and rt-pcr products confirmed the organization of these two genes. both genes have closely-related autosomal or x paralogs although these paralogs are not annotated. homologs of gyg and gyg were not identified in an. stephensi by blastn and tblastx searches of the genome assemblies as well as the trace files, although a short repetitive region of gyg does have significant alignments to the an. stephensi genome from tblastx. a portion of the gyg transcript shares more than percent nucleotide identity with the agap group of genes on both the x and autosomes, while the flanking regions are distantly but significantly related to agap and agap respectively. a portion of the gyg transcript shares greater than percent nucleotide identity with a multi-copy non-coding sequence in the an. gambiae genome. another part of gyg also contains a region that is homologous to agap. the gyg and gyg transcripts are reverse complementary to agap and agap respectively. this transcription orientation coupled with the fact that gyg and gyg are in opposite orientation may indicate functions related to double-stranded rna formation, rather than coding for proteins. because there are no known orthologs to agap agap and agap we could not determine the ancestry of the gyg and gyg genes relative to their paralogs. gyg was discovered in an assembly of our own an. gambiae illumina data. there are two similar sequences in the an. gambiae genome, both sharing percent nucleotide identity with gyg. one is on chromosome 2 l and the other is not mapped but its chromosome quotient suggests it is an autosomal or x sequence. thus, gyg is again involved in a recent duplication, although we cannot assign the direction of the duplication. a blastx search shows that gyg is distantly related to agap, which belongs to a large family of general transcription factor ii-i repeat domain-containing proteins. discussion here we have shown that the cq method can reliably identify y sequences using the ratio of alignments from male and female sequence data. using the cq method we identified y chromosome sequences encompassing more than bases of novel y chromosome sequences from the malaria mosquitoes an. stephensi and an. gambiae. it is also possible to identify y sequences by coverage depth. in male sequence data, y sequences should have half the coverage of autosomal sequences. in mixed male and female sequence data, y sequences should have one-quarter the coverage of autosomal sequences. the cq method does not directly rely on coverage depth to identify y sequences. while y sequences often do have low coverage depth, it is often not sufficient to classify a sequence as y-linked solely by low coverage depth. mis-assemblies of the reference genome, and allelic differences between the reference genome and the short-read sequence data could lead to low coverage depth causing errant classification of y-linkage. unmasked repetitive sequences may lead to high coverage depth, obfuscating many y sequences. instead, the cq method uses the ratio of female to male alignments to a reference sequence. the cq method takes into account both the lack of female alignments and the presence of male alignments when classifying sequences as y-linked, reducing the overall rate of false positives and false negatives. mis-assemblies and allelic variation can be detected by the lack of male reads aligned to a sequence. for a sequence to be classified as y-linked with the chromosome quotient it must have many alignments from males, and zero or a few alignments from females. although the functions of the y genes that we found are not yet clear, there are some tantalizing hints to their potential functions. syg is homologous to an hmg-box transcription factor involved in a number of developmental functions. hmg-box transcription factors include the two key genes involved in male-determination in mammals. syg shares homology with a general transcription factor ii repeat domain containing protein, indicating the possibility that it is a y-linked transcription factor. five of the y genes are clearly expressed in the testes suggesting that they may perform male-specific roles. all six y genes are expressed in the early embryo, which is the time when sexual differentiation begins. our research offers an opportunity to investigate y chromosome evolution, which until now has been restricted to limited lineages. primate y chromosomes have been shown to undergo strict conservation following rapid gene loss. in d. melanogaster, however, almost all the y genes are the result of duplications to the y. here we have shown that two of the three an. stephensi y genes were recently acquired from the autosomes or the x. all three an. gambiae y genes have highly similar paralogs on the autosomes or x, again consistent with recent duplication. thus, gene acquisitions through duplications appear to be a major source of y genes in anopheles mosquitoes. this is intriguing given that unlike d. melanogaster, the anopheles y is thought to harbor a dominant male-determining factor that initiates sexual differentiation, similar to humans and several non-drosophila flies. so far, no conservation between the genes on the two anopheles y chromosomes has been found, although the two species belong to the same subgenus. it is possible that we simply missed the common y genes that originated prior to the divergence between an. gambiae and an. stephensi. however, such shared y genes are less likely to have highly similar autosomal or x paralogs compared to recent duplications, and they may even be unique to the y. y genes that are unique or significantly diverged from their autosomal or x paralogs are the easiest to find with the cq method. to mitigate the effect of potential underrepresentation of y genes in the anopheles genome sequences, we compared the trace sequence files of both species and found no indications of shared y genes. thus, the evidence so far suggests that the an. gambiae and an. stephensi y chromosomes are very different in their gene content. this conclusion is consistent with the observation that two of the three y genes in an. stephensi were recently acquired after the divergence of the two mosquitoes. further investigation of y genes in additional anopheles species will reveal whether the evolutionary scenarios demonstrated by the drosophila y chromosomes, namely recent y gene acquisition and y chromosome replacement, are widely applicable in species with a dominant male-determining factor on the y. some of the y genes described in this study have unique y chromosome sequences. because of their exclusive paternal inheritance and general lack of recombination, these y sequences could serve as powerful molecular markers for the investigation of population structure and incipient speciation, which are often relevant to vectorial capacity. only female mosquitoes bite and transmit pathogens. genetic strategies to control mosquito-borne diseases may include the release of sterile males, the release of insects carrying a dominant lethal gene, or the replacement of a pathogen-susceptible population with a pathogen-resistant one. under these scenarios, the release of males is either required or preferred. a better understanding of the possible functions of mosquito y genes in sexual differentiation will facilitate the production of sterile males for improved sit. it may also enable the production of highly competitive males to improve sit and ridl and to help spread pathogen resistance genes.  in this study, we have shown that y sequences can be identified by their distinctive chromosome quotients, the ratio of female to male alignments. although the genomes of many organisms with y chromosomes have been sequenced, y sequences have been characterized from only a small subset of these genomes. y sequences are surely present in these genomes, but are probably fragmented into short contigs challenging the traditional experimental methods of finding y sequences. the quality of the genome assembly is not a limiting factor in finding y sequences, making the cq method widely applicable to species with fragmented genome. autosome or x sequences rarely have chromosome quotients less than, giving the cq method a low false positive rate. the cq method is able to identify the majority of known y sequences from species we have analyzed. we have also shown that the cq method can identify recent duplications to the y that share more than percent nucleotide identity with autosomal or x paralogs. coupling the cq method with transcriptome sequencing allowed us to identify six y genes in an. stephensi and an. gambiae. analysis of the six y genes indicates rapid y chromosome evolution between an. stephensi and an. gambiae. the anopheles y genes characterized in this study provide unique markers for population and phylogenetic analysis, and opportunities for novel mosquito control measures through the manipulation of sexual dimorphism and fertility. the identification of y genes in more organisms will help us understand the evolution and functions of y chromosomes. 
using the cq method, we successfully identified known y sequences from homo sapiens and drosophila melanogaster. the chromosome quotient is a novel approach to systematically discover y chromosome genes. all six genes are expressed in the early embryo. the anopheles y genes discovered by the cq method provide unique markers for population and phylogenetic analysis, and opportunities for novel mosquito control measures through the manipulation of sexual dimorphism and fertility. y chromosomes are responsible for the initiation of male development, male fertility, and other male-related functions in diverse species. the cq method can reliably identify y chromosome sequences using the ratio of alignments from male and female sequence data. comparisons to transcriptome sequence data with blastn led to the discovery of six anopheles y genes, three from each species. the cq method is widely applicable to species with fragmented genome assemblies produced from next-generation sequencing data. however, y genes are rarely characterized outside a few model species due to the arduous nature of studying the repeat-rich y. stephensi y genes were recently acquired from the autosomes or the x. although an. gambiae belong to the same subgenus, we found no evidence of y genes shared between the species. analysis of the six y genes characterized in this study indicates rapid y chromosome evolution between an. the female to male ratio of the number of alignments to a reference sequence, a parameter called the chromosome quotient, is used to determine whether the sequence is y-linked. in the cq method, genomic dna from males and females is sequenced independently and aligned to candidate reference sequences. the cq method facilitated the discovery of y chromosome sequences from the malaria mosquitoes anopheles stephensi and an.
 y chromosomes are responsible for the initiation of male development, male fertility, and other male-related functions in diverse species. however, y genes are rarely characterized outside a few model species due to the arduous nature of studying the repeat-rich y.  the chromosome quotient is a novel approach to systematically discover y chromosome genes. in the cq method, genomic dna from males and females is sequenced independently and aligned to candidate reference sequences. the female to male ratio of the number of alignments to a reference sequence, a parameter called the chromosome quotient, is used to determine whether the sequence is y-linked. using the cq method, we successfully identified known y sequences from homo sapiens and drosophila melanogaster. the cq method facilitated the discovery of y chromosome sequences from the malaria mosquitoes anopheles stephensi and an. gambiae. comparisons to transcriptome sequence data with blastn led to the discovery of six anopheles y genes, three from each species. all six genes are expressed in the early embryo. two of the three an. stephensi y genes were recently acquired from the autosomes or the x. although an. stephensi and an. gambiae belong to the same subgenus, we found no evidence of y genes shared between the species.  the cq method can reliably identify y chromosome sequences using the ratio of alignments from male and female sequence data. the cq method is widely applicable to species with fragmented genome assemblies produced from next-generation sequencing data. analysis of the six y genes characterized in this study indicates rapid y chromosome evolution between an. stephensi and an. gambiae. the anopheles y genes discovered by the cq method provide unique markers for population and phylogenetic analysis, and opportunities for novel mosquito control measures through the manipulation of sexual dimorphism and fertility. 