This is a test version of Biostars. For the public version, visit https://www.biostars.org.
Can't load reference sequence from file 'GRCh37.fa': Unexpected character 'M' found.

Hi all:

I am trying to do an analysis using GRCh37.fa as reference genome. After running command

pureclip -i aligned.f.duplRm.pooled.R2.bam -bai aligned.f.duplRm.pooled.R2.bam.bai -g GRCh37.fa -iv "1;2;3;4;5;6;7;8;9;10;11;12;13;14;15;16;17;18;19;20;21;22;X;Y;" -nt 10 -o PureCLIP.crosslink_sites.bed

I received an error:

ERROR: Can't load reference sequence from file 'GRCh37.fa': Unexpected character 'M' found.

I got advice from the developer as:

The problem is coming from an external library which is used and which expects the reference sequence to contain only the letters 'A', 'C', 'G', 'T' or 'N'. I know it is not ideal, but if you convert all non-ACGTs to Ns, the problem should be solved

Does anyone can teach me how to convert all non-ACGT to Ns so that I will be able to give it a try?

Thanks,

next-gen-sequencing rna-seq
• 2,623 views • link
updated by Ram Ram • written by Kai_Qi Kai_Qi

It is indeed strange that your reference contains the letter M. As a first step I would double check that the reference contains nucleotides and not amino acids. Once you are sure that this is the case you can use pyfasta of some other tool depending in which programming language you are proficient.

• 0 views • link
written by Fabio Marroni Fabio Marroni

It is a reference genome. Thank you for your comment. I used a differential reference genome and it generated bed file. Nevertheless, I could not see peaks when load the bed to IGV. I guess I need to ask around.

Thanks,

• 0 views • link
written by Kai_Qi Kai_Qi

1 answer

Hi Kai, Got the same error in Pureclip using the reference genome fasta file downloaded from UCSC. To get the clean reference genome, you can try this:

cat hg38as.fa | sed '/^\>/! s/[RYKMSWBVHD]/N/g' > hg38as_clean.fa

It basically replaces any letters in the character set [RYKMSWBVHD] with an N, only in the lines that are not the header aka do not start with >. This fixed the Pureclip error for me.

• 0 views • link
written by swathishekaran93 • swathishekaran93

Thank you, I will try it, but I ended up using ctk toolkit to get it work.

• 0 views • link
written by Kai_Qi Kai_Qi

Log in to answer this question.

More posts like this