Create Test Data for Kiwifruit Genomics Workflows

In [1]:

    

module load samtools/1.3.1


    

In [2]:

    

module load bedtools/2.21.0


    

In [4]:

    

REFFILE=/powerplant/workspace/cfljam/HighHealth/PoolSeq/PS1_reference/PS1_1.68.5.fasta
INDEXFILE=/powerplant/workspace/cfljam/HighHealth/PoolSeq/PS1_reference/PS1_1.68.5.fasta.fai
head $INDEXFILE


    

    




CHR1	19474402	86	19474402	19474403
CHR2	15865742	19474575	15865742	15865743
CHR3	23687156	35340404	23687156	23687157
CHR4	12906829	59027647	12906829	12906830
CHR5	18272032	71934563	18272032	18272033
CHR6	17395956	90206682	17395956	17395957
CHR7	20282415	107602725	20282415	20282416
CHR8	25307543	127885227	25307543	25307544
CHR9	16512530	153192857	16512530	16512531
CHR10	19412578	169705478	19412578	19412579

In [5]:

    

awk -v OFS='\t' '{print $1,1,1000}' $INDEXFILE | head


    

    




CHR1	1	1000
CHR2	1	1000
CHR3	1	1000
CHR4	1	1000
CHR5	1	1000
CHR6	1	1000
CHR7	1	1000
CHR8	1	1000
CHR9	1	1000
CHR10	1	1000

In [6]:

    

bedtools getfasta -fi ${REFFILE} -bed <(awk -v OFS='\t' '{print $1,1,1000}' ${INDEXFILE}) -fo ./kiwitest.fasta


    

In [8]:

    

grep '>' kiwitest.fasta | head


    

    




>CHR1:1-1000
>CHR2:1-1000
>CHR3:1-1000
>CHR4:1-1000
>CHR5:1-1000
>CHR6:1-1000
>CHR7:1-1000
>CHR8:1-1000
>CHR9:1-1000
>CHR10:1-1000

In [9]:

    

wgsim


    

    




Program: wgsim (short read simulator)
Version: 0.3.2
Contact: Heng Li <lh3@sanger.ac.uk>

Usage:   wgsim [options] <in.ref.fa> <out.read1.fq> <out.read2.fq>

Options: -e FLOAT      base error rate [0.020]
         -d INT        outer distance between the two ends [500]
         -s INT        standard deviation [50]
         -N INT        number of read pairs [1000000]
         -1 INT        length of the first read [70]
         -2 INT        length of the second read [70]
         -r FLOAT      rate of mutations [0.0010]
         -R FLOAT      fraction of indels [0.15]
         -X FLOAT      probability an indel is extended [0.30]
         -S INT        seed for random generator [0, use the current time]
         -A FLOAT      discard if the fraction of ambiguous bases higher than FLOAT [0.05]
         -h            haplotype mode

In [14]:

    

wgsim -N 1000 -1 125 -2 125 -r 0.01 -S 10 ./kiwitest.fasta kiwitest.1.R1.fq kiwitest.1.R2.fq > kiwitest1.wgsim.out
wgsim -N 1000 -1 125 -2 125 -r 0.01 -S 5 ./kiwitest.fasta kiwitest.2.R1.fq kiwitest.2.R2.fq > kiwitest2.wgsim.out


    

    




[wgsim] seed = 10
[wgsim_core] calculating the total length of the reference sequence...
[wgsim_core] 30 sequences, total length: 29970
[wgsim] seed = 5
[wgsim_core] calculating the total length of the reference sequence...
[wgsim_core] 30 sequences, total length: 29970

In [15]:

    

head kiwitest1.wgsim.out


    

    




CHR1:1-1000	245	G	T	-
CHR1:1-1000	288	T	W	+
CHR1:1-1000	394	C	M	+
CHR1:1-1000	406	G	-	-
CHR1:1-1000	500	C	-	-
CHR1:1-1000	541	C	T	-
CHR1:1-1000	617	T	G	-
CHR1:1-1000	926	G	C	-
CHR1:1-1000	929	T	G	-
CHR1:1-1000	965	T	-	-

In [16]:

    

rm *.gz
pigz kiwitest*.fq


    

In [17]:

    

ls -lh


    

    




total 665K
-rw-rw-r--. 1 cfphxd cfphxd 8.1K Aug 18 10:59 2016-07-06 Create Test Data.ipynb
-rw-rw-r--. 1 cfphxd cfphxd  36K Aug 18 10:59 kiwitest.1.R1.fq.gz
-rw-rw-r--. 1 cfphxd cfphxd  36K Aug 18 10:59 kiwitest.1.R2.fq.gz
-rw-rw-r--. 1 cfphxd cfphxd 6.8K Aug 18 10:59 kiwitest1.wgsim.out
-rw-rw-r--. 1 cfphxd cfphxd  36K Aug 18 10:59 kiwitest.2.R1.fq.gz
-rw-rw-r--. 1 cfphxd cfphxd  36K Aug 18 10:59 kiwitest.2.R2.fq.gz
-rw-rw-r--. 1 cfphxd cfphxd 7.1K Aug 18 10:59 kiwitest2.wgsim.out
-rw-rw-r--. 1 cfphxd cfphxd  30K Aug 18 10:56 kiwitest.fasta
-rw-rw-r--. 1 cfphxd cfphxd  939 Aug 16 15:00 kiwitest.fasta.fai
-rw-rw-r--. 1 cfphxd cfphxd 6.3K Jul  7 09:40 kiwitest.wgsim.out
-rw-rw-r--. 1 cfphxd cfphxd  146 Jul  7 09:40 Readme.md

In [1]:

    

zcat kiwitest.1.R1.fq.gz | head


    

    




@CHR1:1-1000_471_914_3:1:1_6:0:0_0/1
TCCAACTGGGCCGAGACTTGGGTGACCCGATGCTTTGTCCACATTATCCTCTGATGGATCGACATGCTTTGTCCGTGGTACAGTCTTACCATTGTCTGAACAAGACTCAATCTTCTCGTCGAAGA
+
22222222222222222222222222222222222222222222222222222222222222222222222222222222222222222222222222222222222222222222222222222
@CHR1:1-1000_246_703_4:0:0_3:1:0_1/1
TCCACTCTTCTTTAAATTGGTAATGGCATTAGTAGTTTTGTTTAAAGATGAATAAATAGATGAATCATATTTTCCTATGAACCGACCAGAGGTCGATACAAAGCAAGTCTACGAGGAGAAGATTG
+
22222222222222222222222222222222222222222222222222222222222222222222222222222222222222222222222222222222222222222222222222222
@CHR1:1-1000_273_830_3:2:0_3:0:0_2/1
CAAAAGTAGTCATGGTTGCAACTTGCAATGATATTATCTAAATGCTACTTTTCTTGAAATGTGGAATCTACCATATTGTCTAAGAATCATGGTCTGTACTTTTAGTAATTTAATATTTGAAGAAA

gzip: stdout: Broken pipe

In [2]:

    

zcat kiwitest.2.R1.fq.gz | head


    

    




@CHR1:1-1000_46_657_4:2:0_3:1:0_0/1
AGCTGCACTCTACAGATGAGACCACATTGTTGAAAACTAAAGCGTCACAAGAAATTTCATATTACCACCAATAAAAATAAAAGTAAAAAAAGAACGAAAAAGATAGAACATATAACTTTCCCAAC
+
22222222222222222222222222222222222222222222222222222222222222222222222222222222222222222222222222222222222222222222222222222
@CHR1:1-1000_43_610_4:2:0_2:2:0_1/1
CGATACAAAGCAAGTCTTGGAGGAGAAGATTGAGACTTGTTCAGACAATGGTAATCCTGTACCACGGACGAAGCATGTCGATCCATCAGAGGATAATGTGGACAAAGCATGGGGGTCAGCCAAGT
+
22222222222222222222222222222222222222222222222222222222222222222222222222222222222222222222222222222222222222222222222222222
@CHR1:1-1000_270_744_1:0:0_2:0:0_2/1
GAATGCAATCCACTGTATTTCAATTGTAAACAGTCACAACATTATGGCGGATGGAAAGAGTTCACTTGTTGCTCATGGTGGCTGTTGACTGGTAAGGCTAACCACAGGTACATCCAAGCTTCTCC

gzip: stdout: Broken pipe

In [ ]: