Merging Pairwise Sequence Alignment

[loba17]

Dear All,

I am looking for a way to "merge" two large sequence files. The first sequence is the reference sequence and the second is a consensus sequence containing ambiguous nucleotides (Ns). The merged sequence should correspond to the consensus sequence but the missing nucleotides should be used from the reference.

Example:

BAM-consensus: CCCTANNNNNNNNCATCTACATGG
Reference: GTATAGATATCATCATCTACATCC
Output => CCCTAGATATCATCATCTACATGG


I tried Emboss (cons, consambig, megamerger) but the problems are ambiguous nucleotides. I wrote a simple unix script and it worked just fine but not with large sequences (100MB). Any ideas that would work for larger files?

Thanks for considering my question!

[RickBioinf]

Maybe try and right a Perl file, I think that would be your best option for this.
I personally don't know a tool for this.

[gsgs]

m1:a=fgetc(file1);
b=fgetc(file2);
if(a=='N')a=b;
fputc(a,file3);
if(feof(file1)==0)goto m1;

[loba17]

thanks gsgs

for your helpful reply. unfortunately I am not familiar with fgetc. I am on a Linux RedHat system.

[gsgs]

it's a pity with all those different systems ...
do you have a C-compiler ? The simpler the better

[loba17]

the system was not my choice

I do have a c compiler an it is :
gcc version 4.4.6 20120305 (Red Hat 4.4.6-4)

does this help in anyway?

[gsgs]

yes, that could work in theory ...

have you ever successfully compiled a simple C-program ?

I'll test it (on Windows) and post it here later (or PM or email if you prefer)

are the 2 files in fasta-format ?

[loba17]

I guess I would be able to handle compiling a simple C program using:
gcc infile.c -o outfile

I have one file with two sequences but I can easily split into two.

I could also try to compile the lines you send me myself if you could help me with the C formating:

/* File: merger.c
Function: Simple program to merge a reference sequence and with another sequence.
*/

#include <stdio.h>

void main ()
{
m1:a=fgetc(file1);
b=fgetc(file2);
if(a=='N')a=b;
fputc(a,file3);
if(feof(file1)==0)goto m1;
}

???

Thanks again for your help!

[gsgs]

/* I'm not sure about Linux, (ascii 10 for end of line rather than
ascii 13 and ascii 10 for windows
it prints the result, redirect to a file via ">" - piping :
loba oldfile.fa ref.fa > newfile.fa

fasta-file with one sequence only, same length and lines and aligned to the reference
should it also handle multi-sequence-fasta-files ?

*/


#include <stdio.h>
int a,b;
FILE *file1,*file2;

int main(int argc,char*argv[]){
if(argc<3){printf("\nusage:loba file1 file2\n\n");
printf("replaces N in file1 by what file2 has at that position\n");
exit(1);}


if((file1=fopen(argv[1],"rb"))==NULL){printf("\ncan't open file %s\n",argv[1]);exit(1);}
if((file2=fopen(argv[2],"rb"))==NULL){printf("\ncan't open file %s\n",argv[2]);exit(1);}

a=0;while(a!=10){a=fgetc(file1);printf("%c",a);} /* first line printed */
a=0;while(a!=10)a=fgetc(file2); /* first line ignored (fasta-format ?)*/

m1:;
a=fgetc(file1);
b=fgetc(file2);
if(a=='N')a=b;
printf("%c",a);
if(feof(file1)==0)goto m1;
fclose(file1);fclose(file2);

}

[loba17]

it works - thanks !

the compiler finishes with warning but no errors (see below) - I guess that is ok?
I tested it and it works now I have to test it with a big file !!

Thank you very much for your help !

merger.c: In function ‘main’:
merger.c:21: warning: incompatible implicit declaration of built-in function ‘exit’
merger.c:23: warning: incompatible implicit declaration of built-in function ‘exit’
merger.c:24: warning: incompatible implicit declaration of built-in function ‘exit’

[gsgs]

if it works, the warnings are no problem
maybe it should be exit() instead of exit(1) for your compiler