2025-02-27
El formato fasta está basado en texto.
Se utiliza representar secuencias bien de nucleótidos bien de aminoácidos.
Tanto unos como otros son representados por una sola letra.
También tiene símbolos para representar un hueco (gap) o parada en la traducción o bien que no se sabe el nucleótido o aminoácido.
Tiene una línea que comienza con el símbolo \(>\) al que sigue una descripción de la secuencia.
En la siguiente línea empieza la secuencia de bases o aminoácidos. Se recomienda que no tener más de 80 columnas y se pueden tener todas las filas que se precisen.
Un ejemplo
>MCHU - Calmodulin - Human, rabbit, bovine, rat, and chicken MADQLTEEQIAEFKEAFSLFDKDGDGTITTKELGTVMRSLGQNPTEAELQDMINEVDADGNGTIDF PEFLTMMARKMKDTDSEEEIREAFRVFDKDGNGYISAAELRHVMTNLGEKLTDEEVDEMIREADI*
El formato fastq es el más popular para datos de secuencias.
Consiste de cuatro líneas por lectura.
La primera que comienza con el carácter @ y contiene el nombre de la secuencia con alguna descripción opcional de la misma.
La segunda línea contiene la secuencia con las letras que correspondan dependiendo del tipo (nucleótidos, aminoácidos).
La tercera línea que comienza con \(+\) contiene información opcional sobre la secuencia.
La cuarta línea cuantifica la confianza o calidad en la determinación de cada base recogida en la segunda línea (Phred).
¿Cómo se cuantifica la confianza o precisión?
x = sample(Biostrings::DNA_ALPHABET[1:4],10,replace=TRUE)
y = runif(10,min=0,max=.001)
names(y) = x
y C T T T C A
0.0005638712 0.0005860426 0.0003363001 0.0004159278 0.0005068044 0.0007981135
T C T G
0.0007861088 0.0003369066 0.0001463833 0.0008107285 Q valores.Tomamos una lectura o secuencia corta.
Pretendemos encontrar en una secuencia larga donde es similar, en que punto hay una mayor similitud respecto de la secuencia larga o de referencia.
Alienamiento global o local.
Utilizamos : Bowtie2 o STAR.
SRA Run SelectorSelect podemos seleccionar información de toda la muestra o bien de una parte.Accession List generamos un fichero texto SRR_Acc_List.txt en donde aparecen los nombres de los Run. Lo guardamos, por ejemplo, con dicho nombre.Metadata el fichero SraRunTable.txt. Contiene las variables fenotípicas o covariables para cada una de las muestras. Lo bajamos y guardamos con el mismo nombre.prefetch de SRA-Toolkit.~/ncbi/public/sra/.PRJNA218851/sra.fastq-dump de SRA-Toolkit.donde files es el fichero SRR_Acc_List.txt, es decir, los nombres de los ficheros sin la extensión.
tami::StarSamtools_sh().genomeDir0 por el directorio donde queremos guardar el fichero de índices y GRCh38.p13.genome.fa por el nombre del fichero fastq donde tenemos el genoma de referencia.f_1.fastq y f_2.fastq y name_file es el nombre del fichero de salida.indexDir0 y queremos alinear el fichero name_file.fastq.name_file.name_file por cada nombre original.tami::StarSamtools_sh().pacman::p_load(Rsamtools,GenomicFeatures,GenomicAlignments,org.Hs.eg.db)
gtfFile = "GRCh38.p13.genome/gencode.v37.annotation.gff3"
txdb = makeTxDbFromGFF(gtfFile, format="gff3")
genes = exonsBy(txdb, by="gene")
indir = getwd()
files = list.files(indir, pattern = '*sort.bam')
bamLst = BamFileList(files, index=character(),obeyQname=TRUE)
PRJNA411984 = summarizeOverlaps(features = genes,
read=bamLst,
mode="Union",
singleEnd=FALSE,
ignore.strand=TRUE,
fragments=FALSE)
save(PRJNA411984,file="PRJNA411984.rda")tamidata::PRJNA266927tamidata::PRJNA297664tamidata::PRJNA297798SummarizedExperimenttamidataassay?[1] "matrix" "array" ICR1 LSR1 NME1 RDN5-1 RDN5-6 RNA170
4 5 6 42 47 50 DataFrame with 6 rows and 4 columns
SampleName Run treatment replication
<character> <character> <factor> <numeric>
1 GSM1900735 SRR2549634 Wild 1
2 GSM1900737 SRR2549636 Wild 3
3 GSM1900739 SRR2549638 SEC66 deletion 2
4 GSM1900736 SRR2549635 Wild 2
5 GSM1900738 SRR2549637 SEC66 deletion 1
6 GSM1900740 SRR2549639 SEC66 deletion 3treatment con[1] Wild Wild SEC66 deletion Wild SEC66 deletion
[6] SEC66 deletion
Levels: Wild SEC66 deletionDataFrame with 6 rows and 4 columns
ORF SGD ENTREZID ENSEMBL
<character> <character> <character> <character>
15S_rRNA 15S_rRNA NA NA NA
21S_rRNA 21S_rRNA NA NA NA
HRA1 HRA1 S000119380 9164866 NA
ICR1 ICR1 S000132612 9164906 NA
LSR1 LSR1 S000006478 9164871 NA
NME1 NME1 S000007436 9164967 NAGRangesList object of length 1:
$YBR031W
GRanges object with 1 range and 2 metadata columns:
seqnames ranges strand | exon_id exon_name
<Rle> <IRanges> <Rle> | <integer> <character>
[1] II 300166-301254 + | 224 YBR031W.1
-------
seqinfo: 17 sequences (1 circular) from an unspecified genome; no seqlengths\(Y_{ij}\) (\(y_{ij}\)) el conteo aleatorio (observado) del gen \(i\) en la muestra \(j\).
Denotamos \(L_i\) la longitud del gen \(i\) y \(m_j\) es el total de lecturas de la librería \(j\) (o tamaño de la librería \(j\)).
En (Robinson and Oshlack 2010) se asume que la media de \(Y_{ij}\) verifica \[
E \bigg [Y_{ij} \bigg ] = \frac{\mu_{ij} L_i}{c_j} m_j,
\]
siendo \(c_j = \sum_{i=1}^N \mu_{ij} L_i\). Notemos que \(c_j\) nos está representando el total de RNA en la muestra.
La producción total de RNA, \(c_j\), no es conocida.
Podemos estimar el cociente de estos valores para dos muestras \[ f_j = \frac{c_j}{c_{j'}} \]
Elegimos una muestra como muestra de referencia. Por ejemplo, la muestra \(r\) denota a partir de ahora la muestra tomada (arbitrariamente) como de referencia.
Nos fijamos en la muestra \(j\) y vamos a determinar la constante por la que multiplicaremos los conteos originales.
En lo que sigue tanto \(j\) como \(r\) son fijos y las cantidades definidas dependen de \(i\) que denota el gen.
Se define \[ M_{ij}^{(r)} = \log_2 \frac{y_{ij}/m_j}{y_{ir}/m_r} = \log_2 (y_{ij}/m_j) -\log_2 ( y_{ir}/m_r), \] y \[ A_{ij}^{(r)} = \frac{1}{2} \bigg ( \log_2 (y_{ij}/m_j) + \log_2 (y_{ir}/m_r) \bigg ) \]
Se eliminan los valores extremos tanto de los \(M_{ij}^{(r)}\) como de los \(A_{ij}^{(r)}\).
En concreto eliminamos un porcentaje de los \(M_i\) más pequeños y el mismo porcentaje de los más grandes.
Lo mismo hacemos para los valores \(A_i\).
También eliminamos aquellos índices \(i\) tales que $y_{ij} =0 $ o bien \(y_{ir} =0\).
El conjunto de índices i restante lo denotamos por \(G^*\).
Finalmente, el factor de normalización sería \[ \log_2(TMM_j^{(r)}) = \frac{\sum_{i \in G^*} w_{ij}^{(r)} M_{ij}^{(r)}}{\sum_{i \in G^*} w_{ij}^{(r)}} \] con \[ w_{ij}^{(r)} = \frac{m_j - y_{ij}}{m_jy_{ij}} + \frac{m_r - y_{ir}}{m_r y_{ir}}. \]
La muestra de referencia \(r\) es fija y lo que acabamos de calcular es el factor por el que multiplicamos los conteos originales de la muestra \(j\). Este factor viene dado por \(TMM_j^{(r)}\).
Obviamente tenemos que \(TMM_r^{(r)} = 1\) y esta muestra de referencia no se normaliza.