In [1]:

    

ls ~/Downloads/*.fastq.gz


    

    




ls: cannot access /home/ilya/Downloads/*.fastq.gz: No such file or directory

In [ ]:

    

mv ~/Downloads/*.fastq.gz ../data


    

In [1]:

    

ls -lah ../data


    

    




total 3.5G
drwxrwxr-x 2 ilya ilya 4.0K Dec  9 16:14 .
drwxrwxr-x 6 ilya ilya 4.0K Dec  2 17:47 ..
-rw-rw-r-- 1 ilya ilya  24M Sep 30 08:54 7E2.pileup.gz
-rw-rw-r-- 1 ilya ilya  21M Sep 30 15:49 AAGCTA_R2.pileup.gz
-rw-rw-r-- 1 ilya ilya 779K Nov 18 16:35 ACAGTG.UTR_heatmap.txt.gz
-rw-rw-r-- 1 ilya ilya 765K Nov 18 16:35 ATCACG.UTR_heatmap.txt.gz
-rw-rw-r-- 1 ilya ilya 146M Sep 29 12:04 BJ-HSR1_R1.fastq.gz
-rw-rw-r-- 1 ilya ilya 1.1M Nov 18 16:35 CGATGT.UTR_heatmap.txt.gz
-rw-rw-r-- 1 ilya ilya 500K Sep 16 09:03 contigs.fasta
-rw-rw-r-- 1 ilya ilya 603K Nov  2 15:12 dfm.csv
-rw-r----- 1 ilya ilya  644 Sep 16 09:00 dHSR1.fa
-rw-rw-r-- 1 ilya ilya  92M Dec  2 17:41 ecoli.idx
-rw-rw-r-- 1 ilya ilya 1.4M Dec  9 16:13 Escherichia_coli_str_k_12_substr_mg1655_gca_000801205.GCA_000801205.1.29.cdna.all.fa.gz
-rw-rw-r-- 1 ilya ilya 1.2M Nov 18 16:35 GCCAAT.UTR_heatmap.txt.gz
-rw-r----- 1 ilya ilya 1.4M Oct 14 14:41 GCF_000005845.2_ASM584v2_genomic.fna.gz
-rw-rw-r-- 1 ilya ilya 440K Sep 10 13:08 GCF_000005845.2_ASM584v2_genomic.gff.gz
-rw-rw-r-- 1 ilya ilya 5.4K Sep 23 15:58 gradtimes.txt
-rw-rw-r-- 1 ilya ilya  445 Sep 16 09:00 hHSR-435.fa
-rw-rw-r-- 1 ilya ilya  611 Sep 16 09:00 hHSR.fa
-rw-rw-r-- 1 ilya ilya 4.0K Oct 28 16:04 ROSE1_25.txt
-rw-rw-r-- 1 ilya ilya 4.2K Oct 28 16:04 ROSE1_26.txt
-rw-rw-r-- 1 ilya ilya 4.3K Oct 28 16:04 ROSE1_27.txt
-rw-rw-r-- 1 ilya ilya 4.4K Oct 28 16:04 ROSE1_28.txt
-rw-rw-r-- 1 ilya ilya 4.5K Oct 28 16:04 ROSE1_29.txt
-rw-rw-r-- 1 ilya ilya 4.7K Oct 28 16:04 ROSE1_30.txt
-rw-rw-r-- 1 ilya ilya 5.0K Oct 28 16:04 ROSE1_31.txt
-rw-rw-r-- 1 ilya ilya 5.2K Oct 28 16:04 ROSE1_32.txt
-rw-rw-r-- 1 ilya ilya 5.7K Oct 28 16:04 ROSE1_33.txt
-rw-rw-r-- 1 ilya ilya 6.0K Oct 28 16:04 ROSE1_34.txt
-rw-rw-r-- 1 ilya ilya 6.5K Oct 28 16:04 ROSE1_35.txt
-rw-rw-r-- 1 ilya ilya 6.7K Oct 28 16:04 ROSE1_36.txt
-rw-rw-r-- 1 ilya ilya 7.1K Oct 28 16:04 ROSE1_37.txt
-rw-rw-r-- 1 ilya ilya 7.8K Oct 28 16:04 ROSE1_38.txt
-rw-rw-r-- 1 ilya ilya 8.1K Oct 28 16:04 ROSE1_39.txt
-rw-rw-r-- 1 ilya ilya 8.2K Oct 28 16:04 ROSE1_40.txt
-rw-rw-r-- 1 ilya ilya 8.6K Oct 28 16:04 ROSE1_41.txt
-rw-rw-r-- 1 ilya ilya 8.8K Oct 28 16:04 ROSE1_42.txt
-rw-rw-r-- 1 ilya ilya 9.3K Oct 28 16:04 ROSE1_43.txt
-rw-rw-r-- 1 ilya ilya  10K Oct 28 16:04 ROSE1_44.txt
-rw-rw-r-- 1 ilya ilya  11K Oct 28 16:04 ROSE1_45.txt
-rw-rw-r-- 1 ilya ilya  11K Oct 28 16:04 ROSE1_46.txt
-rw-rw-r-- 1 ilya ilya  126 Sep 16 09:04 rose.fa
-rw-r----- 1 ilya ilya 534M Dec  2 15:06 S1_AGTCAA_L003_R1_001.fastq.gz
-rw-r----- 1 ilya ilya 585M Dec  2 15:16 S2_AGTTCC_L003_R1_001.fastq.gz
-rw-r----- 1 ilya ilya 537M Dec  2 15:11 S3_ATGTCA_L003_R1_001.fastq.gz
-rw-r----- 1 ilya ilya 487M Dec  2 15:43 W1_CGATGT_L003_R1_001.fastq.gz
-rw-r----- 1 ilya ilya 555M Dec  2 15:38 W2_TGACCA_L003_R1_001.fastq.gz
-rw-r----- 1 ilya ilya 541M Dec  2 15:43 W3_ACAGTG_L003_R1_001.fastq.gz

In [2]:

    

gunzip -c ../data/Escherichia_coli_str_k_12_substr_mg1655_gca_000801205.GCA_000801205.1.29.cdna.all.fa.gz | head


    

    




>AIZ89464 cdna:novel chromosome:GCA_000801205.1:Chromosome:1941:2132:-1 gene:EO53_00015 transcript:AIZ89464 description:"membrane protein"
ATGATGACCATCAGCGATATCATTGAAATTATTGTCGTTTGCGCACTGATATTTTTCCCG
CTGGGCTATCTGGCGCGGCACTCTTTGCGACGCATTCGCGACACCTTACGTTTGTTCTTT
GCTAAACCTCGTTATGTTAAACCGGCCGGGACGTTACGCCGCACGGAAAAAGCCAGGGCA
ACCAAAAAATGA
>AIZ89465 cdna:novel chromosome:GCA_000801205.1:Chromosome:3972:4160:1 gene:EO53_00025 transcript:AIZ89465 description:"hypothetical protein"
ATGAATAACAATGAACCAGATACTCTGCCTGATCCCGCGATAGGCTATATCTTCCAGAAT
GATATTGTGGCGTTAAAGCAGGCATTTTCACTGCCTGATATTGATTATGCCGATATTTCC
CAACGCGAACAGTTGGCCGCGGCATTAAAACGCTGGCCGTTGCTGGCAGAGTTTGCGCAA
CAAAAGTAG

gzip: stdout: Broken pipe

In [3]:

    

kallisto index -i ../data/ecoli.idx ../data/Escherichia_coli_str_k_12_substr_mg1655_gca_000801205.GCA_000801205.1.29.cdna.all.fa.gz


    

    




[build] loading fasta file ../data/Escherichia_coli_str_k_12_substr_mg1655_gca_000801205.GCA_000801205.1.29.cdna.all.fa.gz
[build] k-mer length: 31
[build] counting k-mers ... done.
[build] building target de Bruijn graph ...  done 
[build] creating equivalence classes ...  done
[build] target de Bruijn graph has 5027 contigs and contains 3851033 k-mers 

In [4]:

    

kallisto quant


    

    




kallisto 0.42.4
Computes equivalence classes for reads and quantifies abundances

Usage: kallisto quant [arguments] FASTQ-files

Required arguments:
-i, --index=STRING            Filename for the kallisto index to be used for
                              quantification
-o, --output-dir=STRING       Directory to write output to

Optional arguments:
    --bias                    Perform sequence based bias correction
-b, --bootstrap-samples=INT   Number of bootstrap samples (default: 0)
    --seed=INT                Seed for the bootstrap sampling (default: 42)
    --plaintext               Output plaintext instead of HDF5
    --single                  Quantify single-end reads
-l, --fragment-length=DOUBLE  Estimated average fragment length
-s, --sd=DOUBLE               Estimated standard deviation of fragment length
                              (default: value is estimated from the input data)
-t, --threads=INT             Number of threads to use (default: 1)
    --pseudobam               Output pseudoalignments in SAM format to stdout

In [5]:

    

ls -lah ../data/*.fastq.gz


    

    




-rw-rw-r-- 1 ilya ilya 146M Sep 29 12:04 ../data/BJ-HSR1_R1.fastq.gz
-rw-r----- 1 ilya ilya 534M Dec  2 15:06 ../data/S1_AGTCAA_L003_R1_001.fastq.gz
-rw-r----- 1 ilya ilya 585M Dec  2 15:16 ../data/S2_AGTTCC_L003_R1_001.fastq.gz
-rw-r----- 1 ilya ilya 537M Dec  2 15:11 ../data/S3_ATGTCA_L003_R1_001.fastq.gz
-rw-r----- 1 ilya ilya 487M Dec  2 15:43 ../data/W1_CGATGT_L003_R1_001.fastq.gz
-rw-r----- 1 ilya ilya 555M Dec  2 15:38 ../data/W2_TGACCA_L003_R1_001.fastq.gz
-rw-r----- 1 ilya ilya 541M Dec  2 15:43 ../data/W3_ACAGTG_L003_R1_001.fastq.gz

In [6]:

    

# This needs to be run from the file to preserve output

wt="W1 W2 W3"
fasta_dir=../data/

for sample in $wt
do
    echo "Running kallisto quant for sample "$sample
    output=../kallisto/$sample
    if [ ! -d "$output" ]; then
        mkdir -p $output
    fi
    read1=$fasta_dir$(ls $fasta_dir | grep $sample)
    kallisto quant --single -l 200 -s 30 -i ../data/ecoli.idx -o $output -b 100 <(gunzip -c $read1)
done


    

    




Running kallisto quant for sample W1

[quant] fragment length distribution is truncated gaussian with mean = 200, sd = 30
[index] k-mer length: 31
[index] number of targets: 4,391
[index] number of k-mers: 3,851,033
[index] number of equivalence classes: 4,578
[quant] running in single-end mode
[quant] will process file 1: /dev/fd/63
[quant] finding pseudoalignments for the reads ... done
[quant] processed 14,312,155 reads, 11,014,220 reads pseudoaligned
[   em] quantifying the abundances ... done
[   em] the Expectation-Maximization algorithm ran for 55 rounds


Running kallisto quant for sample W2

[quant] fragment length distribution is truncated gaussian with mean = 200, sd = 30
[index] k-mer length: 31
[index] number of targets: 4,391
[index] number of k-mers: 3,851,033
[index] number of equivalence classes: 4,578
[quant] running in single-end mode
[quant] will process file 1: /dev/fd/63
[quant] finding pseudoalignments for the reads ... done
[quant] processed 16,255,960 reads, 12,922,297 reads pseudoaligned
[   em] quantifying the abundances ... done
[   em] the Expectation-Maximization algorithm ran for 190 rounds


Running kallisto quant for sample W3

[quant] fragment length distribution is truncated gaussian with mean = 200, sd = 30
[index] k-mer length: 31
[index] number of targets: 4,391
[index] number of k-mers: 3,851,033
[index] number of equivalence classes: 4,578
[quant] running in single-end mode
[quant] will process file 1: /dev/fd/63
[quant] finding pseudoalignments for the reads ... done
[quant] processed 15,874,534 reads, 12,802,660 reads pseudoaligned
[   em] quantifying the abundances ... done
[   em] the Expectation-Maximization algorithm ran for 60 rounds

In [7]:

    

mut="S1 S2 S3"
fasta_dir=../data/

for sample in $mut
do
    echo "Running kallisto quant for sample "$sample
    output=../kallisto/$sample
    if [ ! -d "$output" ]; then
        mkdir -p $output
    fi
    read1=$fasta_dir$(ls $fasta_dir | grep $sample)
    kallisto quant --single -l 200 -s 30 -i ../data/ecoli.idx -o $output -b 100 <(gunzip -c $read1)
done


    

    




Running kallisto quant for sample S1

[quant] fragment length distribution is truncated gaussian with mean = 200, sd = 30
[index] k-mer length: 31
[index] number of targets: 4,391
[index] number of k-mers: 3,851,033
[index] number of equivalence classes: 4,578
[quant] running in single-end mode
[quant] will process file 1: /dev/fd/63
[quant] finding pseudoalignments for the reads ... done
[quant] processed 15,624,516 reads, 12,582,765 reads pseudoaligned
[   em] quantifying the abundances ... done
[   em] the Expectation-Maximization algorithm ran for 125 rounds


Running kallisto quant for sample S2

[quant] fragment length distribution is truncated gaussian with mean = 200, sd = 30
[index] k-mer length: 31
[index] number of targets: 4,391
[index] number of k-mers: 3,851,033
[index] number of equivalence classes: 4,578
[quant] running in single-end mode
[quant] will process file 1: /dev/fd/63
[quant] finding pseudoalignments for the reads ... done
[quant] processed 17,117,029 reads, 13,695,055 reads pseudoaligned
[   em] quantifying the abundances ... done
[   em] the Expectation-Maximization algorithm ran for 55 rounds


Running kallisto quant for sample S3

[quant] fragment length distribution is truncated gaussian with mean = 200, sd = 30
[index] k-mer length: 31
[index] number of targets: 4,391
[index] number of k-mers: 3,851,033
[index] number of equivalence classes: 4,578
[quant] running in single-end mode
[quant] will process file 1: /dev/fd/63
[quant] finding pseudoalignments for the reads ... done
[quant] processed 15,783,453 reads, 11,651,599 reads pseudoaligned
[   em] quantifying the abundances ... done
[   em] the Expectation-Maximization algorithm ran for 54 rounds

In [8]:

    

ls -lah ../kallisto/W1


    

    




total 1.6M
drwxrwxr-x 2 ilya ilya 4.0K Dec  2 17:50 .
drwxrwxr-x 8 ilya ilya 4.0K Dec  2 17:54 ..
-rw-rw-r-- 1 ilya ilya 1.4M Dec  9 16:20 abundance.h5
-rw-rw-r-- 1 ilya ilya 130K Dec  9 16:18 abundance.tsv
-rw-rw-r-- 1 ilya ilya  273 Dec  9 16:18 run_info.json

In [9]:

    

head ../kallisto/W1/abundance.tsv


    

    




target_id	length	eff_length	est_counts	tpm
AIZ89464	192	21.9082	170	193.016
AIZ89465	189	20.9989	193	228.618
AIZ89466	729	530	77	3.61381
AIZ89467	2619	2420	95	0.97647
AIZ89468	2340	2141	128	1.48711
AIZ89469	1107	908	42	1.15057
AIZ89470	3474	3275	290	2.20261
AIZ89471	1956	1757	360	5.09661
AIZ89472	1287	1088	7098	162.277

In [ ]: