This is a test version of Biostars. For the public version, visit https://www.biostars.org.
How to identify Gap location

To submit assembled genome scaffolds for DDBJ, I have to indicate the regions of sequencing gaps.

for example

scaffold1_cov134        assembly_gap    1647..1712
                        assembly_gap    9101..9259
scaffold3_cov149        assembly_gap    1173..1187

It seems to be very time consuming. Does anyone know how to identify Gap (n) location each scaffold ?

gap scaffold location
• 1,459 views • link
written by Takuma • Takuma

1 answer

Try seqkit locate

$ seqkit locate -P -i -r -G -p 'n+' genome.fa \
    | sed 1d | awk '{print $1"\tassembly_gap\t"$5".."$6}'
scaffold1_cov134        assembly_gap    7..9
scaffold3_cov149        assembly_gap    7..9
scaffold3_cov149        assembly_gap    16..16
scaffold3_cov149        assembly_gap    19..19

Removing duplicated sequence IDs.

$ seqkit locate -P -i -r -G -p 'n+' genome.fa \
    | sed 1d | awk '{print $1"\tassembly_gap\t"$5".."$6}' \
    | awk '{cnt[$1]++; if(cnt[$1]>1)$1=""; print}'
scaffold1_cov134        assembly_gap    7..9
scaffold3_cov149        assembly_gap    7..9
 assembly_gap 16..16
 assembly_gap 19..19
• 0 views • link
written by shenwei356 shenwei356

Thanks, shenwei356. I did it. I appreciate your helpful comment and sepkit!

• 0 views • link
written by Takuma • Takuma

Log in to answer this question.

More posts like this