Convergent genomic signatures of domestication in sheep and goats

Sampling

Domestic sheep (O. aries) and goats (C. hircus) were sampled in Iran (IROA and IRCH groups, respectively) and Morocco (MOOA and MOCH groups, respectively) for a total of 20 animals per group (Supplementary Fig. 6). Estas amostras foram recolhidas entre Janeiro de 2008 e Março de 2012 na parte norte de Marrocos e entre Agosto de 2011 e Julho de 2012 no Noroeste do Irão, no âmbito do projecto europeu Nextgen (Grant Agreement no. 244356), de acordo com os regulamentos éticos da Directiva 86/609/CEE da União Europeia. Clipes auriculares foram recolhidos da parte distal da orelha de animais escolhidos aleatoriamente e imediatamente armazenados em etanol a 96% durante um dia antes de serem transferidos em contas de gel de sílica até à extracção do ADN.

As espécies selvagens muflão asiático (O. orientalis) e Bezoar ibex (C. aegagrus) foram amostradas no Noroeste do Irão dentro do berço de domesticação21,22. Treze muflões asiáticos e 18 tecidos de Bezoar ibex (respectivamente, grupos IROO e IRCA, Suplemento Fig. 6) foram coletados de animais em cativeiro ou recentemente caçados, e de amostras congeladas disponíveis no Departamento de Meio Ambiente do Irã. Esta abordagem de amostragem individual foi concebida para minimizar o potencial de enviesamento, evitando a sobre-representação dos efeitos locais (por exemplo wpOA incluiu 20 amostras de re-sequenciamento de genoma inteiro (WGS) a 12x de cobertura, representando 20 raças mundiais diferentes fornecidas pelo International Sheep Genome Consortium. wpCH consistiu de 14 amostras WGS sequenciadas a 12x de cobertura, representando 9 indivíduos europeus, ou seja 2 amostras de alpinos franceses, e 2 amostras de Saanen franceses sequenciadas pelo INRA, 5 amostras de Saanen italianos fornecidas pela Parco Tecnologico Padano, e 5 indivíduos australianos, ou seja, 2 Boer, 2 Rangeland, e 1 Cashmere amostras fornecidas pelo CSIRO (Dados Suplementares 5).

Produção de dados WGS

>

DNA genômico foi extraído com sucesso de todas as amostras de tecido usando o kit de tecido Macherey Nagel NucleoSpin 96, adaptando o protocolo do fabricante. A amostragem do tecido foi realizada em blocos quadrados MN para obter fragmentos de 25 mg por amostra. Foram preparados três blocos e meio de MN quadrado 96, e a extração foi realizada utilizando um manipulador de Tecan Freedom Evo Liquid, seguindo o protocolo do fabricante. Uma etapa de pré-lise foi realizada para homogeneizar amostras com 180 µl de tampão T1 e 25 µl de proteinase K durante a noite a 56 °C. Para ajustar as condições de ligação, foram adicionados 200 µl de tampão BQ1 e a placa de amostra foi incubada 1 h a 70 °C; posteriormente foram adicionados 200 µl de etanol a 100%. Os lisados foram transferidos para a placa de ligação do tecido Nucleospin Tissue e um vácuo (-0,2 bar, 5 min) foi aplicado para remover o fluxo. Três etapas de lavagem foram feitas com tampões BW e B5, respectivamente, e um vácuo foi aplicado novamente para descartar o fluxo. Antes da eluição do DNA genômico, uma membrana de sílica da placa de ligação do tecido Nucleospin Tissue foi seca sob vácuo com pelo menos -0,6 bar durante 10 min. A etapa de eluição foi realizada com 100 µl de tampão BE pré-aquecido (70 °C) e uma etapa de centrifugação a 3700 rcf durante 5 min em poços de 96-PCR. O DNA genômico foi armazenado a 4 °C para evitar congelamento-degelo e testado para concentração (como ng/μl) usando o método Picogreen e usando um Nanodrop.

Todos os genomas foram ressequenciados de 500 ng de DNA genômico que foram tosquiados a uma faixa de 150-700 bp usando o instrumento Covaris® E210 para cada amostra e usado para preparação da biblioteca Illumina® por um protocolo semi-automatizado. O reparo final, a adaptação A-tailing e os adaptadores compatíveis com Illumina® (BioScientific) foram realizados utilizando o SPRIWorks Library Preparation System e o instrumento SPRI TE (Beckmann Coulter) seguindo o protocolo do fabricante. Foi aplicada uma selecção de tamanho 300-600 bp para recuperar a maior parte dos fragmentos. Os fragmentos de ADN foram amplificados por 12 ciclos de PCR usando o Platinum Pfx Taq Polymerase Kit (Life® Technologies) e os iniciadores específicos do adaptador Illumina®. As bibliotecas foram purificadas com contas 0,8x AMPure XP (Beckmann Coulter) e analisadas com o Bioanalisador Agilent 2100 (Agilent® Technologies) e quantificação qPCR. Bibliotecas foram sequenciadas usando 100 leituras de leitura química de base em célula de fluxo paired-end no Illumina® HiSeq2000.

Illumina paired-end reads for Ovis foram mapeadas para o genoma de referência das ovelhas (OAR v3.1, conjunto GenBank GCA_000298735.146), e para Capra para o genoma de referência das cabras (CHIR v1.0, conjunto GenBank GCA_000317765.147), usando BWA-MEM48. O arquivo BAM produzido para cada indivíduo foi classificado usando Picard SortSam e melhorado usando sequencialmente Picard MarkDuplicates (http://picard.sourceforge.net), GATK RealignerTargetCreator e GATK IndelRealigner49, e SAMtools calmd50.

Variant discovery foi realizado usando três algoritmos diferentes: Samtools mpileup50, GATK UnifiedGenotyper51, e Freebayes52. Os locais das variantes foram identificados independentemente para cada um dos seis grupos, usando os modos multiamostra dos algoritmos de chamada: (i) 162 amostras de MOOA; (ii) 20 amostras de IROA; (iii) 14 amostras de IROO; (iv) 162 amostras de MOCH; (v) 20 amostras de IRCH; (vi) 19 amostras de IRCA. Para alguns grupos, os WGS de mais indivíduos estavam disponíveis como parte do projeto NextGen (ver acima). As amostras utilizadas no presente estudo foram selecionadas para obter grupos equilibrados de 20 indivíduos sempre que possível. Para os grupos IRCA e IROO, amostras adicionais ficaram disponíveis numa fase posterior e foram adicionadas para análises a jusante. Animais com baixo alinhamento e qualidade de chamada foram removidos para obter o conjunto de dados final (Dados Suplementares 5).

Em cada grupo, houve duas rodadas sucessivas de filtragem de qualidade do local variante. O estágio 1 de filtragem fundiu as chamadas dos três algoritmos, enquanto filtrava as chamadas de menor confiança. Um local variante passou se foi chamado por pelo menos dois algoritmos diferentes de chamada com qualidade de variante de fred >30. Um alelo alternativo em um local passou se foi chamado por qualquer um dos algoritmos de chamada, e a contagem do genótipo foi >0. A filtragem do estágio 2 usou a Recalibração de Pontuação de Qualidade de Variante pelo GATK. Primeiro, geramos um conjunto de treinamento dos locais de maior confiança de variantes dentro do grupo onde (i) o local é chamado pelos três chamadores de variantes com qualidade de variante em escala de feixe >100, (ii) o local é bialélico, (iii) a contagem de alelos menores é de pelo menos 3, enquanto a contagem de amostras com qualidade de escala de feixe do genótipo >30. O conjunto de treinamento foi usado para construir um modelo gaussiano usando a ferramenta GATK VariantRecalibrator usando as seguintes anotações de variantes do UnifiedGenotyper: QD, HaplotypeScore, MQRankSum, ReadPosRankSum, FS, DP, InbreedingCoefficient. Um modelo Gaussiano foi aplicado ao conjunto de dados completo, gerando um VQSLOD (log odds ratio de ser uma variante verdadeira). Os sites foram filtrados se VQSLOD <cutoff value. O valor de corte foi definido para cada grupo pelo seguinte: VQSLOD mínimo = {o valor mediano da VQSLOD para variantes do conjunto de treinamento}-3 × {o desvio mediano absoluto da VQSLOD de variantes do conjunto de treinamento}. O rácio SNP de transição/transversão sugeriu que o critério de corte escolhido deu o melhor equilíbrio entre selectividade e sensibilidade.

SNPs foram gerados conjuntos de chamadas para seis grupos de animais Ovis e Capra (i.e., domésticos iranianos e marroquinos, e selvagens para cada género). Como as análises realizadas neste estudo exigiram comparações entre grupos, criamos conjuntos de chamadas genótipos em um conjunto consistente de locais de SNP para todos os animais de qualquer grupo. Para cada gênero, fundimos os locais de chamada de variantes de seus três grupos, e apenas mantivemos posições bialleicas sem dados ausentes. Os genótipos foram chamados novamente em cada local SNP bialélico para todos os indivíduos de interesse pelo GATK UnifiedGenotyper, utilizando a opção GENOTYPE_GIVEN_ALLELES. Nesta fase, a lista de indivíduos foi ampliada para incluir os animais pertencentes aos painéis mundiais de raças ovina e caprina (wpOA e wpCH) e amostras selvagens adicionais que ficaram disponíveis nesta fase (4 O. orientalis e 4 C. aegagrus). Os genótipos foram melhorados e faseados dentro dos grupos pelo Beagle 453, e depois filtrados para fora onde a probabilidade do genótipo era inferior a 0,95. Finalmente, filtramos locais que eram monomórficos nos diferentes subconjuntos de indivíduos utilizados neste estudo (ver abaixo).

Para comparar os sinais de seleção detectados entre Ovis e Capra, realizamos um alinhamento cruzado entre os dois genomas de referência. Primeiro, utilizamos o pipeline de alinhamento par a par do código base54 do Ensembl release 69 para alinhar os genomas de referência de ovinos (OARv3.1) e caprinos (CHIR1.0). Este pipeline usa LastZ55 para alinhar ao nível do DNA, seguido de pós-processamento no qual blocos alinhados são encadeados de acordo com a sua localização em ambos os genomas. O pipeline de alinhamento em pares LastZ é executado rotineiramente pelo Ensembl para todas as espécies suportadas, mas o bode ainda não está incluído no Ensembl. Para evitar o enviesamento para qualquer das espécies, produzimos dois alinhamentos interespecíficos diferentes. Um usava ovelhas como genoma de referência e cabras como não referência, enquanto o outro usava cabras como genoma de referência e ovelhas como não referência. A diferença é que as regiões genómicas das espécies de referência são forçadas a mapear exclusivamente a loci único da espécie não-referencial, enquanto que as regiões genómicas não-referenciais são autorizadas a mapear para múltiplas localizações da espécie de referência. Obtivemos para segmentos de cromossomos de um genoma de referência as coordenadas no genoma não-referencial. Finalmente, para os SNPs descobertos em um gênero, usamos todo o alinhamento do genoma com o genoma de referência do outro gênero para identificar as posições correspondentes (Tabela Complementar 6).

Estrutura genética

A fim de descrever a diversidade genética dentro dos grupos, usamos VCFtools56 para calcular a estatística resumida da variação genética dos 73 indivíduos para Ovis (ou seja 13 IROO, 20 IROA, 20 MOOA e 20 wpOA) e 72 indivíduos para Capra (ou seja, 18 IRCA, 20 IRCH, 20 MOCH, e 14 wpCH). As estatísticas medidas foram o número total de variantes polimórficas (S) para todo o conjunto de indivíduos em cada gênero e dentro de cada grupo, a diversidade média de nucleotídeos (π) dentro de cada grupo e o coeficiente de consanguinidade (F) para cada indivíduo. Dentro de cada gênero, as diferenças entre o grupo selvagem e cada grupo doméstico foram testadas usando um teste t unilateral para valores individuais de consanguinidade e carga genética, e um teste Mann-Whitney de dois lados para diversidade de nucleotídeos por local.

A divergência geral entre os quatro grupos dentro de cada gênero (ou seja selvagens, iranianos e marroquinos domésticos, e painel mundial) foi estimado usando todos os SNPs bialleicos e o par médio ponderado Fst seguindo Weir e Cockerham57 como implementado no VCFtools56. A estrutura genética entre os grupos foi avaliada com o método de agrupamento sNMF26, após a poda do conjunto de dados para remover SNPs com disequilíbrio de ligação (r²) maior que 0,2 usando VCFtools. O desequilíbrio de ligação (r²) foi calculado entre pares de SNPs dentro de janelas deslizantes de 50 SNPs, com um SNP por par removido aleatoriamente quando r² era maior que 0,2. Para cada análise sNMF, cinco séries do mesmo número de clusters (K) foram realizadas com valores de K de 1 a 10. Usamos o critério de centralidade cruzada para identificar a solução de agrupamento mais provável, entretanto, partições alternativas para diferentes números de K também foram exploradas para avaliar como os indivíduos foram divididos entre agrupamentos.

Para separar entre ancestralidade compartilhada e mistura, corremos o TreeMix27 para estimar conjuntamente as divisões da população e subsequentes eventos de mistura usando o conjunto de dados podados usado para sNMF. Corremos o TreeMix com a opção -global para refinar as nossas inferências de máxima verosimilhança. Enraizamos a árvore TreeMix com a divisão entre indivíduos selvagens e domésticos. O tamanho do bloco para o jackknifing foi -k 500 SNPs, o que corresponde aproximadamente a 150 kb, excedendo a média dos blocos de LD encontrados tanto em ovinos como em caprinos. Geramos uma árvore de Máxima Probabilidade sem migração e depois adicionamos eventos de migração e examinamos a mudança incremental na variância explicada pelo modelo e os valores residuais entre os indivíduos. O objetivo foi detectar qualquer potencial alto valor residual ou margem de migração entre os indivíduos selvagens e domésticos. Para explorar melhor a relevância estatística dos possíveis vetores de mistura identificados pela TreeMix (Tabela Suplementar 3), calculamos o teste de três populações f328 como um teste formal de introgressão genética, usando o programa qp3Pop do conjunto ADMIXTOOLS58 para cada combinação de grupos. Para Capra, o grupo wpCH foi dividido entre raças australianas, raças francesas e raças italianas. Os resultados são relatados em Supplementary Data 2.

Demographic inference

Para cada gênero, realizamos análises de inferência demográfica ancestral usando o modelo MSMC implementado no software MSMC225. O MSMC é baseado no coalescente Markoviano sequencialmente em pares59; entretanto, ele usa haplótipos de dados de sequência de genoma em fases como entrada. Para cada análise usamos dois indivíduos de um grupo, portanto 4 haplótipos. Cada análise foi repetida para outro conjunto aleatório de dois indivíduos, ou seja, uma réplica da análise por grupo. Os arquivos de entrada e saída foram gerados e analisados com os scripts python fornecidos com o software MSMC e encontrados em https://github.com/stschiff/msmc-tools. Os parâmetros de análise foram mantidos como padrão, exceto a taxa de mutação que foi definida para 2,5×10-8 e a duração da geração foi definida para 2 anos. Para estimar a incerteza nas estimativas de tempo, variamos esses parâmetros (taxa de mutação de 2,5×10-8 e 1,0×10-8 em combinação com a duração da geração de 2 e 4 anos) e fornecemos uma estimativa aproximada do período de domesticação (ver Figura Complementar 2).

Carga genética

Carga genética foi estimada de duas maneiras. Primeiro, calculando a carga genética de cada indivíduo como a soma dos efeitos deletérios sobre todas as posições genômicas que codificam a proteína, seguindo o método de Librado et al.60. Resumidamente, como um proxy para restrição evolutiva, usamos os escores PhyloP do alinhamento dos mamíferos de 46 vias (http://hgdownload.cse.ucsc.edu/goldenPath/hg19/phyloP46way/placentalMammals/). A partir deste alinhamento, identificamos locais codificadores de proteínas evoluindo sob restrições funcionais (escore de PhyloP ≥1.5). Para cada genoma Ovis ou Capra, investigamos então se esses locais sofreram mutações. Se sim, somamos os escores de phyloP em todos os locais mutantes, de modo que mutações em locais altamente restritos contribuem proporcionalmente mais para a estimativa de carga total. Isto forneceu uma estimativa de carga para cada genoma ovino/caprino. Finalmente, para obter uma carga média por site, dividimos pelo número total de posições analisadas. Vale ressaltar que condicionamos em locais homozigotos para evitar a modelagem do coeficiente de dominância das mutações em locais heterozigotos (por exemplo, recessivos, intermediários, dominantes). Em segundo lugar, comparamos gene por gene a carga deletéria genética em grupos Ovis selvagens e domesticados realizando um teste Wilcoxon, sendo a hipótese alternativa que os animais domésticos têm mais carga que os parentes selvagens. os valores de p foram corrigidos para testes múltiplos61 e aplicamos um limiar de valores de p ajustados < 0,05. Realizamos uma análise de enriquecimento ontológico de genes no conjunto de genes mostrando um aumento significativo na carga genética usando WebGestalt62,63. Como os genomas de referência são pouco anotados para os genes, confiamos em ortologs de cópia única entre nossa espécie e a humana e o rato. Os genes do cromossoma X foram excluídos do conjunto de fundo. Não realizamos esta análise em Capra devido à maior consanguinidade observada nas amostras selvagens.

Detecção de assinaturas de selecção

Para detectar assinaturas de selecção relacionadas com a domesticação, utilizámos todos os SNPs bialleicos mostrando uma frequência alélica menor que 0,10 em pelo menos um dos três grupos testados (i.e., grupos domésticos Iranianos e Marroquinos, e o grupo selvagem para cada género). Como esperávamos que as assinaturas de selecção relacionadas com o processo de domesticação estivessem presentes em todos os animais domésticos, adoptámos a seguinte estratégia geral: testámos com hapFLK29 (ver Nota Complementar 5 e Figuras Complementares 9, 10 e 11) para cada género o grupo selvagem contra cada um dos grupos domésticos tradicionalmente geridos (i.e., iranianos e marroquinos) e concentrámo-nos naquelas regiões comuns sob selecção que foram detectadas em ambos os casos. Os tamanhos das amostras de grupo (n = 13-20) foram compatíveis com os requisitos do método29. Verificamos visualmente se as assinaturas consistentes de seleção encontradas com hapFLK também estavam presentes no conjunto correspondente do painel mundial de cada gênero, mas não incluímos esses grupos no teste estatístico devido à sua composição multi-sangue. Finalmente, procuramos sinais compartilhados de seleção entre Ovis e Capra usando uma abordagem estratificada de FDR. A estratégia está descrita na Fig. 4.

Realizámos testes hapFLK para contrastar o grupo selvagem com cada um dos grupos iranianos e marroquinos de cada género. A matriz de parentesco foi calculada a partir das distâncias genéticas do Reynold64 entre pares de grupos, usando um subconjunto aleatório de um por cento das variantes. A árvore populacional inferida foi construída utilizando o algoritmo de união de vizinhos. Para cada SNP, foi realizado o teste hapFLK que incorpora informações haplotípicas para aumentar o poder de detecção de varreduras seletivas. Para cada SNP testado, a estatística hapFLK calculou o desvio das frequências haplotípicas em relação ao modelo neutro estimado pela matriz de parentesco65. Para explorar a informação de desequilíbrio de ligação, o hapFLK utiliza o modelo multiponto Scheet e Stephens’66 para genótipos multilocus que podem ser ajustados a dados não faseados. Uma das principais aplicações deste modelo é realizar estimativa de fase (software fastPHASE66). Na nossa análise, o modelo foi treinado em dados não faseados e, portanto, a nossa análise contabiliza a incerteza de fase. O método foi usado para reagrupar haplótipos locais ao longo dos cromossomos em um número especificado de clusters K definidos como 25, usando um Modelo Markov Oculto.

Para identificar as regiões comuns que estão putativamente sob seleção nos dois grupos domésticos tradicionalmente gerenciados para cada gênero, combinamos as duas análises hapFLK anteriores. Para cada análise as pontuações hapFLK foram ajustadas a uma distribuição χ2 para obter os valores p (script disponível em https://forge-dga.jouy.inra.fr/projects/hapflk/documents). Os resultados dos dois contrastes entre o grupo selvagem e cada um dos grupos domésticos foram combinados usando o método de Stouffer67 para obter valores de p únicos para a comparação de animais selvagens vs. animais domésticos. Finalmente, a estrutura FDR68 foi aplicada a todo o conjunto de SNPs para converter os p-valores combinados em q-valores. SNPs mostrando valores q < 10-2 foram retidos e agrupados em regiões genômicas quando estavam a menos de 50 kb distantes umas das outras.

Para investigar se o sinal de seleção foi compartilhado entre Ovis e Capra, usamos primeiro o alinhamento cruzado dos dois genomas de referência para identificar segmentos homólogos. Em seguida, aplicamos uma estrutura FDR estratificada69. Esta abordagem baseia-se no fato de que existe uma estratificação inerente nos testes dada a informação prévia nos dados genéticos69, pois a distribuição subjacente das verdadeiras hipóteses alternativas pode ser diferente de acordo com as diferentes dinâmicas das várias regiões genômicas, levando a diferentes distribuições de p-valores. Isto requer a obtenção de valores p ajustados de FDR (isto é, valores q) separadamente para os diferentes estratos. Buscamos convergências em cada gênero, separando as regiões homólogas às detectadas no outro gênero (referido como o estrato compartilhado) e o resto do genoma (referido como o estrato geral). Extraímos os valores p separadamente para cada um dos dois estratos definidos e depois calculamos os valores q através da estrutura FDR. Esses valores q estratificados foram as quantidades finais consideradas para significância estatística (<10-2) para detectar SNPs em seleção e fundi-los nas regiões genômicas correspondentes.

Para testar a convergência de assinaturas de seleção diferenciando animais selvagens de animais domésticos em ambos os gêneros, examinamos a relação entre o limiar de significância aplicado aos valores q (que fizemos variar de 0,2 a 0,002) em um gênero e a probabilidade estimada de que um SNP seja selecionado no estrato compartilhado do outro gênero usando Storey et al.70 abordagem. Um aumento na probabilidade inferida com uma diminuição do limiar aplicado ao valor q (aumento na stringency) indica que quanto mais significativa a região estiver em um gênero, maior a probabilidade de encontrarmos SNPs significativos no outro gênero.

Filtramos os sinais de seleção que não eram consistentes entre os três grupos domésticos. Para cada região detectada, utilizamos os haplótipos faseados de indivíduos que foram agrupados usando árvores Neighbor-Joining com base na porcentagem de identidade entre sequências. Apenas as regiões que apresentaram sinais consistentes foram mantidas (Figura Complementar 5).

Para inferir se os sinais de seleção detectados com hapFLK indicavam relaxamento da seleção ou seleção positiva nos domésticos, estimamos a diferença na diversidade de nucleotídeos (π) em cada região putativa sob seleção entre os grupos selvagens e domésticos. Expressamos esta diferença como o índice Δπ, que foi calculado para cada região genômica como a diferença entre π calculada para o grupo selvagem e a média de π para os grupos domésticos Iranianos e Marroquinos, menos a diferença em π entre estes dois grupos calculada sobre todo o genoma:

$$\Delta \pi = \ esquerda( {\pi _{{\rm wilds}}} – PI _{{{\i1}mathrm{\i}operatorname{\i1}iran-morocco}}}}}} Certo)_{{{\i1}{\i1}-direito)_{\i1}-operatorname{\i1}genomic-region}}}} – esquerda… – PI _{{{\i1}mathrm{\i}operatorname{\i1}iran-morocco}}}}} Certo)_{{{\i1}{\i1}mathrm{\i1}operatorname{\i1}whole-genome}}}}$$

Um valor negativo indicaria que a diversidade de nucleotídeos é menor no grupo selvagem em comparação com a média dos dois grupos domésticos, e seria considerado como mostrando um relaxamento da seleção nestes últimos grupos, diversificando a seleção nos domésticos ou a seleção positiva nos selvagens. Ao contrário, um valor positivo indicaria uma seleção direcional positiva ou estabilizadora que ocorreu nos grupos domésticos. Também utilizamos o agrupamento de haplótipos para verificar manualmente em cada região se a varredura seletiva detectada confirmava as indicações dadas pelo índice Δπ.

Fizemos as interpretações funcionais a seguir. Para cada região em seleção, consideramos a região mais 50 kb em cada lado para identificar papéis funcionais e 5 kb a montante e a jusante dos genes e avaliamos a sobreposição entre essas coordenadas para reter os genes de interesse. Finalmente consideramos que um gene estava relacionado a uma determinada região detectada quando as posições da região e do gene estavam sobrepostas. Avaliamos então qual gene era o mais provavelmente visado pela seleção, considerando o gene mais próximo do sinal superior, ou seja, a posição do valor q mais baixo dentro da região. Os genes foram anotados funcionalmente usando Uniprot (http://www.uniprot.org/), considerando seu envolvimento em 30 termos infantis (ou seja, descendentes diretos dos termos) da categoria “Processo Biológico” (ou seja, GO:0008150). Recuperamos todos os termos GO correspondentes a cada gene (Dados Suplementares 4) para 30 das 33 categorias, pois não consideramos três termos que não estavam envolvidos nas funções dos mamíferos (i.e., GO:0006791 utilização do enxofre, GO:0006794 utilização do fósforo, GO:0015976 utilização do carbono). Realizamos dois testes χ2 para comparar as distribuições de genes nas categorias GO, ou seja, (i) genes sob seleção de regiões específicas do gênero versus os de regiões homólogas, e (ii) todos os genes sob seleção versus os 18.689 genes humanos associados aos termos GO em Swiss-Prot. A fim de interpretar as funções dos genes em um contexto pecuário, também recuperamos a informação disponível na literatura sobre seus efeitos fenotípicos.

Finalmente, para encontrar os SNPs dentro das regiões previamente detectadas que eram as mais diferenciadas entre os grupos selvagens e domésticos, utilizamos a estatística FLK. Quanto ao hapFLK, ele representa o desvio das frequências alélicas de marca única em relação ao modelo neutro estimado pela matriz de parentesco65. O mesmo procedimento foi utilizado para adequar os escores das duas análises a uma distribuição χ2 e combinar os valores de p obtidos como foi utilizado para o teste hapFLK. Entretanto, a distribuição não uniforme dos valores de p impediu a aplicação do quadro FDR e selecionamos SNPs dentro das regiões detectadas com hapFLK mostrando p-valores <10-4. Para estes SNPs utilizamos as anotações Variant Effect Predictor (VEP)71 que foram geradas a partir da anotação do genoma OARv3.1 da ovelha Ensembl v74 para Ovis (http://www.ensembl.org/Ovis_aries/Tools/VEP) e da anotação do genoma da cabra CHIR1.0 produzida pelo gasoduto de anotação do genoma eucariótico da NCBI para Capra (https://www.ncbi.nlm.nih.gov/genome/annotation_euk/process/). Os SNPs foram classificados como intergénicos, a montante e a jusante (incluindo os UTRs), e posições intrónicas e exónicas. As diferenças entre as distribuições de SNPs com valores p FLK <10-4 e todos os SNPs utilizados para detectar assinaturas de seleção foram examinados com um teste χ2-test.

Disponibilidade de dados

Sequências e dados de metadados gerados para as 73 amostras Ovis e 72 Capra utilizadas nestas análises estão disponíveis ao público. Informações gerais e todos os arquivos vcf podem ser encontrados no site da Ensembl (http://projects.ensembl.org/nextgen/). Todos os arquivos Fastq, Bam e de novo assemblages de O. orientalis e C. aegagrus podem ser encontrados no arquivo europeu de nucleotídeos (https://www.ebi.ac.uk/ena) sob o código de acesso do projeto Nextgen (PRJEB7436).

Deixe uma resposta

O seu endereço de email não será publicado.