Initializing a 4-D array

Hello to all,

I need help initializing a 4-D array that will hold all instances from [0,0,0,0] to [4,4,4,4] (only using 1-4). I need it for a comparison chart. What I am trying to do is scan a 5 million character string for each particular 4-mer that occurs and keep count. Any help would be great. Below is the code I have so far.

Code:

---

char t1, t2, t3, t4 ;
        int index ;
        int tupleCount[4][4][4][4] ;

 for( index = 0 ; index < strlen(seqData) - 3 ; ++index )
            {
                t1 = seqData[index] ;
                t2 = seqData[index + 1] ;
                t3 = seqData[index + 2] ;
                t4 = seqData[index + 3] ;
                // count the different 4-mers
            }

---

Sorry, my 4-D array should hold (1,1,1,1) to (4,4,4,4).

Quote:

---

from [0,0,0,0] to [4,4,4,4] (only using 1-4).

---

I do not get it
0-4 is 5 values
using 1-4 values you cannot get [0,0,0,0] vector

so what do you want to do?
how seqData is related to what you need?

to fill 4D array - use 4 loops one inside the other

Code:

---

for(i=0;i<4;i++)
  for(j=0;j<4;j++)
      for(k=0;k<4;k++)
        for(l=0;l<4;l++)
            tupleCount[i][j][k][l] = ...;

---

Quote:

---

Sorry, my 4-D array should hold (1,1,1,1) to (4,4,4,4).

---

Not (0,0,0,0)? Hmm. In that case, you either waste space or subtract one every time you index the array.

Code:

---

for( index = 0 ; index < strlen(seqData) - 3 ; ++index )

---

Computing strlen() is quite expensive. You probably shouldn't do it in a loop. Consider calculating this value beforehand and storing it in a variable.

Now. Unfortunately, this code is broken. :)

Code:

---

                t1 = seqData[index] ;
                t2 = seqData[index + 1] ;
                t3 = seqData[index + 2] ;
                t4 = seqData[index + 3] ;

---

When you use seqData, you get a char ****. When you use seqData[x], you get a char ***. And so on, until you get a char when using seqData[x][y][z][a].

t1, t2, t3, and t4 are all chars. Therefore, you can't the sort of assignment you have.

I'm not really sure what you're trying to do here. Perhaps you could explain what a 4-mer is, for a start.

Sorry for the confusion; I was rushing when I posted. You are right. I should have said (0,0,0,0) to (3,3,3,3).


I guess it will be easier if I elucidate a little. I am trying to scan a sequence of DNA from the first letter to the last letter -3 for all the 4-mers that occur. I need to collect all the 4-mers so that I can see the total distribution. The array is to hold my count of each 4-mer as I progress through the sequence. I was also matching 0->a, 1->c, 2->g, an 3->t. I am not sure if this is a redundant step because when I print my distribution it will contain letters, not numbers. Hopefully this helps.

Of course it does. If you name your array as tupleCount, and use a suggested method to map a,c,g,t to 0,1,2,3, which i will name acgt_to_int() then:

Code:

---

/* TODO: Verify that the seqData pointer points to enough data */
t1 = acgt_to_int(seqData);
t2 = acgt_to_int(seqData + 1);
t3 = acgt_to_int(seqData + 2);
t4 = acgt_to_int(seqData + 3);
/* TODO: Validate that t1..t4 are correctly in range */
tupleCount[t1][t2][t3][t4]++;

---

I leave the TODO tags to you.

Is there a simple c function to translate a,c,g,t to 0,1,2,3 ??

That is up to you to create. Using a function will mask the gory details of your conversion.
You can do something like:

Code:

---

int acgt_to_int(char *x)
{
    if(*x == 'a')
        return 0;
    ....
}

---

Guess the rest. Also note the fact that you should handle case-sensitivity, and find a return value to indicate error.
I'd suggest -1.

I am not sure what I did wrong but, I am getting all kinds of errors and warnings. I've looked through the code and can't pick up on anything. Any help would be great.

Code:

---

// A Program to count the count of 4-mers in a nucleotide sequence.

#include <stdio.h>
#include <string.h>
#include <ctype.h>
#include <stdlib.h>

FILE *input ;
        FILE *outpur ;
        char buffer[1000] ;       
        int i = 0 ;
        int c ;
        int w,x,y,z ;
        char tr ;
       
        char *seqData ;
     
        seqData = (char *)malloc(10000000) ;
       
     
        char t1, t2, t3, t4 ;
        int index ;

int acgt_to_0123(char *tr)
    {
        if(*tr == 'a')
            {
            return 0 ;
            }
        else if(*tr == 'c')
            {
            return 1 ;
            }
        else if(*tr == 'g')
            {
            return 2 ;
            }
        else if(*tr == 't')
            {
            return 3 ;
            }
        else
            {
            printf("Non acgt character\n") ;
            exit(0) ;
            }
    }

main( int argc, char **argv )
    {
               
       
        // Initialize 4-D array
        int tupleCount[4][4][4][4] ;
        for(w = 0 ; w < 4 ; w++ )
            for(x = 0 ; x < 4 ; x++ )
                for(y = 0 ; y < 4 ; y++ )
                    for(z = 0 ; z < 4 ; z++ )
                        {
                            tupleCount[w][x][y][z] ;
                        }
       
        // Open input file to read from       
        if( ! ( input = fopen( argv[1], "r" ) ) )
            {
                printf( "COULD NOT OPEN FILE %s - Exit!\n", argv[1]) ;
                exit(1) ;
            }       
       
       
        // Collect sequence from GenBank file
        while(fgets(buffer, 1000, input))
                        {
                        // start obtaining bases after ORIGIN
                        if(strstr(buffer, "ORIGIN"))
                            {                                                 
                              while((c=getc(input)) != '/' && c != EOF)
                                    {
                                    if(c >= 'a' && c <= 'z')
                                        {
                                        seqData[i++] = c ;
                                        }
                                    }                           
                            }         
                        }
                       
       
       
       
        // Scan DNA sequence for each 4-mer
        for( index = 0 ; index < strlen(seqData) - 3 ; ++index )
            {
                t1 = acgt_to_0123(seqData[index]) ;
                t2 = acgt_to_0123(seqData[index + 1]) ;
                t3 = acgt_to_0123(seqData[index + 2]) ;
                t4 = acgt_to_0123(seqData[index + 3]) ;
               
               
                // Accumulate a count to find distribution
                tupleCount[t1][t2][t3][t4]++ ;
               
            }
   
                 
        fclose(input) ;
   
        printf("Here is the distribution of 4-mers:\n\n%s", tupleCount ) ;   
        free(seqData) ;
        return(0) ;
    }

---

thanks

Code:

---

        char tr ;

---

Unused variable. Remove it.

Code:

---

        seqData = (char *)malloc(10000000) ;

---

It is not allowed. Global variables must be initialized in block scope ( inside main() ).
Why so many globals? Move those declarations to the start of the main() block.

Code:

---

        // Initialize 4-D array
        int tupleCount[4][4][4][4] ;
        for(w = 0 ; w < 4 ; w++ )
            for(x = 0 ; x < 4 ; x++ )
                for(y = 0 ; y < 4 ; y++ )
                    for(z = 0 ; z < 4 ; z++ )
                        {
                            tupleCount[w][x][y][z] ;
                        }

---

Exactly what is that supposed to mean? Forgot the =0 maybe?

Code:

---

        // Scan DNA sequence for each 4-mer
        for( index = 0 ; index < strlen(seqData) - 3 ; ++index )
            {
                t1 = acgt_to_0123(seqData[index]) ;
                t2 = acgt_to_0123(seqData[index + 1]) ;
                t3 = acgt_to_0123(seqData[index + 2]) ;
                t4 = acgt_to_0123(seqData[index + 3]) ;
               
               
                // Accumulate a count to find distribution
                tupleCount[t1][t2][t3][t4]++ ;
               
            }

---

In the definition of your function acgt_to_0123, i proposed a pointer to character to be passed.
You instead pass the character itself. Change the function to work with characters as it seems to be more straightforward to you. (hint: remove all * from acgt_to_0123).

And a sidenote on correctness. Looking at the for loop in your code and read aloud what is going on. You are moving the index, ONE character at a time, but i presume that the 4-mers are 4 characters long? Or is it considered for your excercise valid to do things like:
YOUR INPUT: 'acaactgagatc'
YOUR 4-MERS: 'acaa', 'caac', 'aact', 'actg', 'ctga', 'tgag', 'gaga', 'agat', 'gatc'

Thanks for the input; will implement. What out said for the 4-mer is exactly right. The input sequence is 5 million bases and all of those examples I must collect. For instance:

count of 4-mers:

aaaa = 3
aaac = 2
.
.
.
tttt = 7

thats what i need printed out at the end.

Again, thanks for the help. I really appreciate it. When I am done this semester I will make sure to go back and get more fluent in the language.

Okay, I've gotten it to compile but now it will not stop running. Did I enter an endless loop somewhere?

Exactly how big is the file being scanned? Post bytes. Because the call to strlen on every iteration of the for loop that does statistics, will do say if N is the size of the file N^2 memory accesses. It should finish, but if N is 1E+6, then after 1E+12 memory accesses. And i think N is bigger in your case. Compute the result of strlen(seqData) in a temporary variable. And save yourself the quadratic performance over N to linear. What an optimization !;)

Code:

---

int lengthDNA = strlen(seqData) ;
       
        // Scan DNA sequence for each 4-mer
        for( index = 0 ; index < (lengthDNA - 3) ; ++index )
            {
                t1 = acgt_to_0123(seqData[index]) ;
                t2 = acgt_to_0123(seqData[index + 1]) ;
                t3 = acgt_to_0123(seqData[index + 2]) ;
                t4 = acgt_to_0123(seqData[index + 3]) ;
               
               
                // Accumulate a count to find distribution
                tupleCount[t1][t2][t3][t4]++ ;               
            }

---

The file I used as input was 11MB and the length of the sequence is 5 Million chars. Still keeps running endlessly

Do a verification printf on lengthDNA, before the loop.

Code:

---

int lengthDNA = strlen(seqData) ;
       
        printf("here is the lenght of the DNA: %d", lengthDNA) ;
       
        // Scan DNA sequence for each 4-mer
        for( index = 0 ; index < (lengthDNA - 3) ; ++index )
            {
                t1 = acgt_to_0123(seqData[index]) ;
                t2 = acgt_to_0123(seqData[index + 1]) ;
                t3 = acgt_to_0123(seqData[index + 2]) ;
                t4 = acgt_to_0123(seqData[index + 3]) ;
               
               
                // Accumulate a count to find distribution
                tupleCount[t1][t2][t3][t4]++ ;               
            }

---

It stopped the loop and printed out the sentences. Why was that able to stop the loop?
Also, it did not print out the distribution of 4-mers.

Code:

---

printf("Here is the distribution of 4-mers:\n\n%s", tupleCount ) ;

---

Quote:

---

Originally Posted by gkoenig

The file I used as input was 11MB and the length of the sequence is 5 Million chars. Still keeps running endlessly

---

Then you may want to think about reading file in big chunks instead of byte-byte

and also avoid using function like strlen on arrays of 5M chars

Is something like this in the ballpark of the initialization thingy?

Code:

---

#include <stdio.h>

int main(void)
{
  static const char dna[] = "acgt";
  char pattern[4 * 4 * 4 * 4][5];
  int w,x,y,z, i = 0;
  for ( w = 0; w < 4; ++w )
  {
      for ( x = 0; x < 4; ++x )
      {
        for ( y = 0; y < 4; ++y )
        {
            for ( z = 0; z < 4; ++z )
            {
              pattern[i][0] = dna[w];
              pattern[i][1] = dna[x];
              pattern[i][2] = dna[y];
              pattern[i][3] = dna[z];
              pattern[i][4] = '\0';
              printf("&#37;3d %s\n", i, pattern[i]);
              ++i;
            }
        }
      }
  }
}

---

I checked the variable seqData and it was able to print out all 5 million bases before.

Here is my code:

Code:

---

// A Program to count the count of 4-mers in a nucleotide sequence.

#include <stdio.h>
#include <string.h>
#include <ctype.h>
#include <stdlib.h>


int acgt_to_0123(char x)
    {
        if(x == 'a')
            {
            return 0 ;
            }
        else if(x == 'c')
            {
            return 1 ;
            }
        else if(x == 'g')
            {
            return 2 ;
            }
        else if(x == 't')
            {
            return 3 ;
            }
        else
            {
            printf("Non acgt character\n") ;
            exit(0) ;
            }
    }

main( int argc, char **argv )
    {
           
        FILE *input ;
        FILE *outpur ;
        char buffer[1000] ;       
        int i = 0 ;
        int c ;
        int w,x,y,z ;
       
        char *seqData ;
     
        seqData = (char *)malloc(10000000) ;
       
     
        char t1, t2, t3, t4 ;
        int index ;
 
       
        // Initialize 4-D array
        int tupleCount[4][4][4][4] ;
        for(w = 0 ; w < 4 ; w++ )
            for(x = 0 ; x < 4 ; x++ )
                for(y = 0 ; y < 4 ; y++ )
                    for(z = 0 ; z < 4 ; z++ )
                        {
                            tupleCount[w][x][y][z] = 0 ;
                        }
       
        // Open input file to read from       
        if( ! ( input = fopen( argv[1], "r" ) ) )
            {
                printf( "COULD NOT OPEN FILE %s - Exit!\n", argv[1]) ;
                exit(1) ;
            }       
       
       
        // Collect sequence from GenBank file
        while(fgets(buffer, 1000, input))
                        {
                        // start obtaining bases after ORIGIN
                        if(strstr(buffer, "ORIGIN"))
                            {                                                 
                              while((c=getc(input)) != '/' && c != EOF)
                                    {
                                    if(c >= 'a' && c <= 'z')
                                        {
                                        seqData[i++] = c ;
                                        }
                                    }                           
                            }         
                        }
                       
       
        int lengthDNA = strlen(seqData) ;
       
        printf("here is the lenght of the DNA: %d", lengthDNA) ;
       
        // Scan DNA sequence for each 4-mer
        for( index = 0 ; index < (lengthDNA - 3) ; ++index )
            {
                t1 = acgt_to_0123(seqData[index]) ;
                t2 = acgt_to_0123(seqData[index + 1]) ;
                t3 = acgt_to_0123(seqData[index + 2]) ;
                t4 = acgt_to_0123(seqData[index + 3]) ;
               
               
                // Accumulate a count to find distribution
                tupleCount[t1][t2][t3][t4]++ ;               
            }
   
                 
        fclose(input) ;
   
        printf("Here is the distribution of 4-mers:\n\n%s", tupleCount ) ;   
        free(seqData) ;
        return(0) ;
    }

---

Code:

---

        while(fgets(buffer, 1000, input))
                        {
                        // start obtaining bases after ORIGIN
                        if(strstr(buffer, "ORIGIN"))
                            {                                                 
                              while((c=getc(input)) != '/' && c != EOF)
                                    {
                                    if(c >= 'a' && c <= 'z')
                                        {
                                        seqData[i++] = c ;
                                        }
                                    }                           
                            }         
                        }

---

1) Nice reading of data INTO buffer 1000 characters from the file.
2) Ok, you find the ORIGIN in buffer and discard it although you certainly NEED it, because the actual data is AFTER the string ORIGIN.
3) getc? Meaning you read 1000 characters and then forget them? Shouldn't you operate on buffer to read the characters from into seqData?

And also the seqData should be 0'ed or else strlen, will be WRONG when you are done. I am not sure how you got after the IO part, in my computer it dumps core right there.
Fix those issues and i am almost certain you are done.

Sorry for the confusion. I am looking through a file for the tag "LOCUS" using strstr and fgets.

Once i hit that tag, I begin collecting all DNA bases with getc directly from the input file. I stop collecting when I hit '/'

I am not having any problem collecting the bases. When I check the contents of seqData, it contains all the bases; all 5 million.

I think my issue is after this point, when I am trying to scan it for all occurrences of 4-mers.

Not sure if this helps but when I print out my tupleCount i get "-1073745432"

Note: When i try to read my characters from buffer instead of input I get the warning:

"warning: passing argument 1 of ‘getc’ from incompatible pointer type"

Thanks

then you may to write your reading loop in less obfuscated way:

Code:

---

        while(fgets(buffer, 1000, input))
        {
                // start obtaining bases after ORIGIN
                if(strstr(buffer, "ORIGIN"))
                { 
                        break; /* end of file header */
                }         
        }
        while((c=getc(input)) != '/' && c != EOF)
        {
                if(c >= 'a' && c <= 'z')
                {
                        seqData[i++] = c ;
                }
        }

---

close file just after you stopped the reading
do not use strlen - use i that specifies how many elements of the array were initialized

%s is not a suatable format for printitng 4D array - use loops and %c format

But vart the code you posted ignores the fact that there are actually data in the buffer, after the header tag. Shouldn't they be copied to the start of seqData? The copy statements should be placed between the two loops.

How do I printout each occurrence in my multidimensional array.

For instance:

0000 = 3
0001 = 2
.
.
.
3333 = 7

The way the code is now I get: "-1073745444" as output.

Here is the code:

Code:

---

// Initialize 4-D array
        int tupleCount[4][4][4][4] ;
        for(w = 0 ; w < 4 ; w++ )
            for(x = 0 ; x < 4 ; x++ )
                for(y = 0 ; y < 4 ; y++ )
                    for(z = 0 ; z < 4 ; z++ )
                        {
                            tupleCount[w][x][y][z] ;
                        }

 
        // Scan DNA sequence for each 4-mer
        for( index = 0 ; index < (lengthDNA - 3) ; ++index )
            {
                t1 = acgt_to_0123(seqData[index]) ;
                t2 = acgt_to_0123(seqData[index + 1]) ;
                t3 = acgt_to_0123(seqData[index + 2]) ;
                t4 = acgt_to_0123(seqData[index + 3]) ;
               
               
                // Accumulate a count to find distribution
                ++tupleCount[t1][t2][t3][t4] ;               
            }
     
        printf("Here is the distribution of 4-mers:\n\n%d\n\n", tupleCount ) ;

---

Thanks for the help

Think about how you zeroed every element of the multidimensional array. The solution is right in front of you ;)

Sorry, but I do not see the answer. Do I have to use a for loop to print out all the occurrences?
If so, how do I define the my range.

for (i = 0 ; i < ????? ; i++)
{
printf("%d",tupleCount[i]
}

????

Hmmm, your frequencies are stored in a 4-dimensional array as you declare.
In your code, ok maybe it's not your code, you use 4-for loops to count over every dimension and set the element to 0.
The idea is the same.

Remodel your problem to printing out all the words from 'aaaa' to 'zzzz'. How would you do it?

Think you have a chessboard and need to visit every square. Think you have a cube made of chessboards and have to visit every square, and finally think you have a hypercube of cubes made of chessboards and have to visit every square, to get confused a bit.

Welcome to the 4th dimension!

All the above problems are equivalent to your problem.

Thanks for the help. Everything works fine. I appreciate that you do not write out the code for me; its important that I ween myself from these forums.

I am maybe too late for this suject. I hope that there is someone who can reply to my question.

I just want to know if i create dynamically a 4D array with pointer like this
char ****pt.
How can I initialize it? can we use the same way with the statical declaration?

there is no difference from 2D arrays

search forum, this question was already asked several times