Começo com um arquivo bam classificado e indexado ("mapped.bam") que representa o mapeamento de pequenas leituras em um genoma de referência e um arquivo de base ("genes.bed") contendo as coordenadas de um conjunto de características de interesse (digamos que sejam genes), para o qual desejo calcular um perfil médio usando programas de deeptools. Eu gostaria de entender as etapas envolvidas para ter certeza do que o eixo vertical do perfil final representa.

Primeira etapa: fazer um arquivo bigwig

Eu crio um arquivo bigwig ("mapped.bw") a partir do arquivo bam usando bamCoverage como segue:

  bamCoverage -b mapped.bam -bs 10 -of = bigwig -o mapped.bw  

A ajuda de bamCoverage diz:

A cobertura é calculada como o número de leituras por bin, onde os bins são janelas curtas de contagem consecutiva de um tamanho definido.

No meu caso, os bins têm 10 bp de comprimento. Minhas leituras são mais longas do que isso.

Para um determinado compartimento, uma determinada leitura pode:

1. sobrepor completamente o compartimento

2. sobrepor o compartimento em n bp, n < 10

3. não sobrepor o compartimento de forma alguma

Por favor, corrija-me se eu estiver errado: Meu palpite é que a leitura é contada como 1 nos casos 1. e 2., e 0 caso contrário, e eu também suponho que uma leitura pode ser contada por vários bins sucessivos se for longa o suficiente .

Segunda etapa: cálculo da média dos genes e plotagem

Eu calculo uma "matriz de meta perfil" ("mapped_on_genes.gz") usando regiões da escala computeMatrix da seguinte maneira:

  computeMatrix scale-regiões \ -S mapped.bw \ -R genes.bed \ --upstream 300 \ --unscaled5prime 500 \ --regionBodyLength 2000 \ --unscaled3prime 500 \ --downstream 300 \ -out mapped_on_genes.gz  

(Existe um parâmetro -bs cujo valor padrão é 10 de acordo com a ajuda do comando.)

Eu uso isso para traçar um perfil usi ng plotProfile :

  plotProfile -m mapped_on_genes.gz \ -out mapped_on_genes_meta_profile.pdf  

Eu obtenho um perfil com valores no eixo y. Em quais unidades estão esses valores?

Meu palpite é o seguinte:

Para o upstream (300 bp) e o 5-primo interno (500 bp), já que o tamanho do compartimento era o mesmo em bamCoverage e computeMatrix , cada ponto no eixo x provavelmente representa uma janela de 10 bp, e sua coordenada y é a média sobre as regiões presentes no arquivo de cama do bins correspondentes no arquivo bigwig, portanto, é um número médio de leituras sobrepondo um bin de 10 bp.

A mesma coisa no lado 3-primo e no lado inferior.

Para os 100 centrais Porção bp, antes de calcular a média sobre as regiões, algum encolhimento ou espalhamento dos compartimentos deve ter sido executado, eu acho que fazendo a média entre os compartimentos vizinhos. Portanto, a unidade final ainda é um número de leituras sobrepondo um compartimento de 10 bp .

E se eu usar compartimentos maiores, devo terminar com valores proporcionalmente mais altos.

Estou correto?