# AULA 3 - ANÁLISE DE TRANSCRIPTOMA
## Segunda 14/10
#### MANHÃ
## Primeira coisa remove a pasta q eu tinhacriado errado rm – r exlista 1
Agora nano exlista1.ssh .
#!/bin/bash
Agora cria o diretório mkdir exlista 1
Cd lista 1
Ai da o wget e cola o link
Ls – lh para ver o tamanho do arquvivo os detalhes
gunzip -c ./*.gz para descompactar
gunzip -c ./*.gz > /exlista1/gene_info.fa
* Primeira coisa qdo abrir o terminal n esquecer de digitar nano e o nome da pasta q quero
Assim eu vou para outra pagina q fica salva, e qdo quiser rodar outra coisa é só alterar nela.
#!/bin/bash #pra fazer o script
mkdir exlista1 #criando o diretório
cd exlista1 #trocando o diretório para esse agr
wget ftp://ftp.ncbi.nih.gov/gene/DATA/gene_info.gz #para baixar um arquivo da net
ls - lh gene_info.gz #para ver o tamanho do arquivo
Para descompactar
gunzip -c gene_info.gz > gene_info #esse mantem a identificação original , mantem os dois arquivos
gunzip gene_info.gz #o problema é q ele n mantem a identificação orginal
grep - P '^7227\t'gene_info > Dm-gene_info.txt
grep – v o que eu não quero
Esse chapéu no começo ^ mostra que eu quero os que tenham 7227 no inicio. Mas isso n basta pq posso trer alguma identificação com 7227 e algum numero dps, então preciso estabeler um limite, por isso se usa o /t
P maiúsculo é pra que essa sintase com o /t funcione. AO invés de usar o – P posso usar o caracter tab, fazendo control v + tab
mkdir – p info/
mv Dm-gene_info.txt info #move esse arquivo para o info
mkdir - p analise/
echo "para trocar de diretorio"
cd analise/
echo "criando link simbolico atalho caminho absoluto"
ln - s /home/gfrezarim/exlista1/info/Dm-gene_info.txt
cut -f 2,3,9 info/Dm-gene_info.txt > Dm-geneID2Symbol.txt #seleciona corta as colunas 2, 3 e 9 do arquivo dm gene que está em info e redireciona pra dm gene id 2 symbol omo pedido no exercício
ftp://hammer.fcav.unesp.br/
login: transcriptomics
senha Un35P
winscp: para passar arquivos do servidor para computador
Bit.ly/bioinfotranscriptomics
Eutilis NCBI quick . Entres programming
Vai no ncbi, vai em gene, e e scolhe um gene de um organismo . Para filtrar o resultado pode clicar em advanced que fica logo abaixo.
All fieds são 3 retangulos, no primeiro coloca gene name, o nome do gene, segundo muda p organismo e coloca o organismo o terceiro n preenche.
O gene name é mt restritivo melhor n usar. Pega o gene, na pagina vou procurar pelo transcrito, na barra la embaixo onde ta escrito mrna and protein. Nm é pra gene codificador, clica no numero de acesso, confere q é o gene e que é mrna.
Ai no primeiro comando eu substituo o numero de acesso que eu quero e vou redirecionar pra transcriptoma . fa . O primeiro fica p sobestrcrever então fica >, os próximos coloca dois >> pra acrescentar
Hspb1 NM_001025569.1
https://www.ncbi.nlm.nih.gov/gene/516099
capn1 NM_174259.2
https://www.ncbi.nlm.nih.gov/nuccore/NM_174259.2
Calpastatina: NM_001030318.3
https://www.ncbi.nlm.nih.gov/gene/281039
DNAJA1 NM_001015637.1
https://www.ncbi.nlm.nih.gov/gene/528862 #acabei ncolocando
IGF-1 NM_001077828.1
https://www.ncbi.nlm.nih.gov/gene/281239
NEDD4 XM_003586572.5
https://www.ncbi.nlm.nih.gov/gene/507781
Primeiro roda isso no diretório msm com login gfrezarim senha fc@v...
> cd ~
> /bin/bash
> perl -MNet::FTP -e \
> '$ftp = new Net::FTP("ftp.ncbi.nlm.nih.gov", Passive => 1);
> $ftp->login; $ftp->binary;
> $ftp->get("/entrez/entrezdirect/edirect.tar.gz");'
> gunzip -c edirect.tar.gz | tar xf -
> rm edirect.tar.gz
> builtin exit
> export PATH=${PATH}:$HOME/edirect >& /dev/null || setenv PATH "${PATH}:$HOME/edirect"
> ./edirect/setup.sh
#isso aqui em cima foi p baixar o programa que tava dando err oqdo a gt tava rodando antes
Isso em cima foi no terminal , agora digita nano pratica1 e coloca isso aquiembaixo.
esearch -db nucleotide -query NM_174259.2 | efetch \
-format fasta >> transcriptoma.fa
esearch -db nucleotide -query NM_001030318.3 | efetch \
-format fasta >> transcriptoma.fa
esearch -db nucleotide -query NM_001015637.1 | efetch \
-format fasta >> transcriptoma.fa
esearch -db nucleotide -query NM_001077828.1 | efetch \
-format fasta >> transcriptoma.fa
esearch -db nucleotide -query XM_003586572.5 | efetch \
-format fasta >> transcriptoma.fa
Antes de rodar precisa habilitar o script pra isso tem q fazer chmod a+x
E rodou, pra saber se deu certo, precisa usar a função grep para ver o que tem na pasta
Então faz grep ‘^>’ transcriptoma . fa e confere q está tudo dentro da pasta
O grep faz pra conferir se a sequencia está la dentro. Eu estipulei que tem q começar com maior > pq toda sequencia começa com maior.
Tem uma função for que facilita o que foi feito, ao invés de fazer um por mim, fz uma função só
É esperado que a cada linha seja executado o acc
rm -f transcriptoma.fa
for acc in NM_057444.3 NM_057439.2 NM_078778.5 NM_057963.5 NM_001032098.2 NM_001032096.2 ; do
echo "Pegando FASTA para ${acc} ..."
esearch -db nucleotide -query ${acc} | efetch \
-format fasta >> transcriptoma.fa
done
Editando para o meu
rm -f transcriptoma.fa
for acc in NM_174259.2 NM_001030318.3 NM_001015637.1 NM_001032098.2 NM_001032096.2 ; do
echo "Pegando FASTA para ${acc} ..."
esearch -db nucleotide -query ${acc} | efetch \
-format fasta >> transcriptoma.fa
done
Foi feito uma pasta chamado sim . Mkdir sim
Então movi o transcriptoma.fa para a pasta sim mv transcriptoma.fa sim/
Cd sim para mudar o diretório para sim
$ grep '^>' transcriptoma.fa
>NM_174259.2 Bos taurus calpain 1 (CAPN1), mRNA
>NM_001030318.3 Bos taurus calpastatin (CAST), transcript variant I, mRNA
>NM_001015637.1 Bos taurus DnaJ heat shock protein family (Hsp40) member A1 (DNAJA1), mRNA
>NM_001077828.1 Bos taurus insulin like growth factor 1 (IGF1), mRNA
>XM_003586572.5 PREDICTED: Bos taurus neural precursor cell expressed, developmentally down-regulated 4, E3 ubiquitin protein ligase (NEDD4), transcript variant X2, mRNA
[gfrezarim@hammer ~]$ nano pratica1.sh
[gfrezarim@hammer ~]$ mkdir sim
[gfrezarim@hammer ~]$ ls
[gfrezarim@hammer ~]$ mv transcriptoma.fa sim/
[gfrezarim@hammer ~]$ ls
E da um ls pra ver se deu certo mover.
Nano abundanceA e nano abundanceB, não esquecer de dar control X pra sair e Y para salvar.
Criando dois arquivos abundance a e abundance b , a coluna a é o numero de acesso o msm que foi usado p criar o trasncriptoma fa, e após o tab na coluna a coloca a proporção .
São as proporções que eu crio no abundance a e b, qual a proporção que mostra de cada sequencia. Pra uma vou mostrar uma certa quantidade, p outra menos. Vamos supor q tem diferença de expressão da condição A para condição B, por ex são simulados duas condições. Alguns genes se mantem na msm prorpoção. Tem q deixar um gene mais expresso em A, do que em B, e tem q deixar uma situação q n tem diferença de expressão.
[gfrezarim@hammer ~]$ cd sim
[gfrezarim@hammer sim]$ ls
transcriptoma.fa
[gfrezarim@hammer sim]$ ^C
[gfrezarim@hammer sim]$ nano abundanceA.txt
[gfrezarim@hammer sim]$ nano abundanceB.txt
[gfrezarim@hammer sim]$ mv pratica1.sh sim
Vao ser gerados 25mil sequencias para cada replica. Esse eh o processo de amostragem do fragmento que será sequenciamento. Vamos recuperar fragmentos de transcritos, com uma média de 300pb. Os fragmentos são gerados por processos mecânicos ou enzimáticos e n se gera fragmentos com exatamente esse tamanho, então existe um certo desvio tolerável, q é o desvio de 30pb.
Renomeia como sample a1 qdo for referente a amostra replica 1
cat ./tmp.frags.r1.1.fasta | renameSeqs.pl \
-if FASTA \
-of FASTA \
-p SAMPLEA1 \
-w 1000 | \
sed 's/^>\(\S\+\).*/>\1/' \ #arquivo fasta sempre começa com > entao n muda
> ./frags_A1.fa
#sed editor de stream, usa pra fazer substituição, tudo o que está nele, por ex no arquivo , quero remover tudo oque n é identificador. O padrão e que se inicia com padrão de maior, tudo que n é espaço mais primeiro lado da expressão é o que quero procurar, o segundo lado é q quero substituir que é tudo que é do sinal de maior, coloca barra e sinal de parênteses, tudo o q ta dentro é capturado.
O sequenciamento pode ser de forma single end ou paired end. Nos dois processos tenho cDNA fragmentados. O sequenciamento single end vai sequenciar apenas um lado, já o paired tanto extremidade 5’3’ e 3’5’ , então o fragmento de 300pb terei então um seqienciamento com 150 de um lado e o outro com mais 150 pb
cat ./frags_A1.fa | simNGS -a \
AGATCGGAAGAGCACACGTCTGAACTCCAGTCACATCACGATCTCGTATGCCGTCTTCTGCTTG:AGATCGGAAGAGCGTCGTGTAGGGAAAGAGTGTAGATCTCGGTGGTCGCCGTATCATT \ #sequencia do adaptador
-p paired \
/usr/local/bioinfo/simNGS/data/s_4_0099.runfile \
-n 151 > ./SAMPLEA1.fastq 2> SAMPLEA1.err.txt #2> capturando erros
Qdo recebemos a sequencia ela vem intercalada. Por ex são duas leituras pq são paired end
Copio o item de 1 a 6, eu faço ABundance 1 e 2, B1 e B2; não esquecer de mudar pra todos os A e B
As vezes se não esta consegijndo abrir o script pode estar em diretório errado. Por ex estava no raw, ai coloca cd ..
Cd ~ vai para o meu home
Cd sim para ir p pasta sim dps q eu tiver no meu home
generate_fragments.py -r transcriptoma.fa \
-a ./abundanceB.txt \
-o ./tmp.frags.r1 \
-t 25000 \
-i 300 \
-s 30
cat ./tmp.frags.r1.1.fasta | renameSeqs.pl \
-if FASTA \
-of FASTA \
-p SAMPLEB2 \
-w 1000 | \
sed 's/^>\(\S\+\).*/>\1/' \
> ./frags_B2.fa
cat ./frags_B2.fa | simNGS -a \
AGATCGGAAGAGCACACGTCTGAACTCCAGTCACATCACGATCTCGTATGCCGTCTTCTGCTTG:AGATCGGAAGAGCG$
-p paired \
/usr/local/bioinfo/simNGS/data/s_4_0099.runfile \
-n 151 > ./SAMPLEB2.fastq 2> SAMPLEB2.err.txt
mkdir -p ./raw
deinterleave_pairs SAMPLEB2.fastq \
-o ./raw/SAMPLEB2_R1.fastq \
./raw/SAMPLEB2_R2.fastq
rm -f ./tmp.frags.r1.1.fasta ./frags_B2.fa ./SAMPLEB2.fastq ./SAMPLEB2.err.txt
wc -l ./raw/SAMPLEB2_R1.fastq
Então eu editei essr arquivo de cima que foi a copia do item de 1 a 6, 4 vezes, uma pra cada abundance que tinha o A e o B então ficou A1 e a2, e b1 e b2 .
Entao tudo isso q foi feito foi pra fragmentar as sequencias em sequencias 300pb, e 25mil sequencias.