code_c/Maasha/src/lib/seq.c

   1 #include "common.h"
   2 #include "seq.h"
   3
   4
   5 void uppercase_seq( char *seq )
   6 {
   7     /* Martin A. Hansen, May 2008 */
   8
   9     /* Uppercase a sequence in place. */
  10
  11     size_t i;
  12
  13     for ( i = 0; seq[ i ]; i++ ) {
  14         seq[ i ] = toupper( seq[ i ] );
  15     }
  16 }
  17
  18
  19 void lowercase_seq( char *seq )
  20 {
  21     /* Martin A. Hansen, May 2008 */
  22
  23     /* Lowercase a sequence in place. */
  24
  25     size_t i;
  26
  27     for ( i = 0; seq[ i ]; i++ ) {
  28         seq[ i ] = tolower( seq[ i ] );
  29     }
  30 }
  31
  32
  33 void revcomp_dna( char *seq )
  34 {
  35     /* Martin A. Hansen, May 2008 */
  36
  37     /* Reverse complement a DNA sequence in place. */
  38
  39     complement_dna( seq );
  40     reverse( seq );
  41 }
  42
  43
  44 void revcomp_rna( char *seq )
  45 {
  46     /* Martin A. Hansen, May 2008 */
  47
  48     /* Reverse complement a RNA sequence in place. */
  49
  50     complement_rna( seq );
  51     reverse( seq );
  52 }
  53
  54
  55 void revcomp_nuc( char *seq )
  56 {
  57     /* Martin A. Hansen, May 2008 */
  58
  59     /* Reverse complements a nucleotide sequence in place. */
  60
  61     complement_nuc( seq );
  62     reverse( seq );
  63 }
  64
  65
  66 void complement_nuc( char *seq )
  67 {
  68     /* Martin A. Hansen, May 2008 */
  69
  70     /* Complements a nucleotide sequence, */
  71     /* after guess the type. */
  72
  73     if ( is_dna( seq ) ) {
  74         complement_dna( seq );
  75     } else if ( is_rna( seq ) ) {
  76         complement_rna( seq );
  77     } else {
  78         die( "Complement nuc failed.\n" );
  79     }
  80 }
  81
  82
  83 void complement_dna( char *seq )
  84 {
  85     /* Martin A. Hansen, May 2008 */
  86
  87     /* Complements a DNA sequence including */
  88     /* ambiguity coded nucleotides. */;
  89
  90     size_t i;
  91
  92     for ( i = 0; seq[ i ]; i++ )
  93     {
  94         switch ( seq[ i ] )
  95         {
  96             case 'a': seq[ i ] = 't'; break;
  97             case 'A': seq[ i ] = 'T'; break;
  98             case 'c': seq[ i ] = 'g'; break;
  99             case 'C': seq[ i ] = 'G'; break;
 100             case 'g': seq[ i ] = 'c'; break;
 101             case 'G': seq[ i ] = 'C'; break;
 102             case 't': seq[ i ] = 'a'; break;
 103             case 'u': seq[ i ] = 'a'; break;
 104             case 'T': seq[ i ] = 'A'; break;
 105             case 'U': seq[ i ] = 'A'; break;
 106             case 'm': seq[ i ] = 'k'; break;
 107             case 'M': seq[ i ] = 'K'; break;
 108             case 'r': seq[ i ] = 'y'; break;
 109             case 'R': seq[ i ] = 'Y'; break;
 110             case 'w': seq[ i ] = 'w'; break;
 111             case 'W': seq[ i ] = 'W'; break;
 112             case 's': seq[ i ] = 'S'; break;
 113             case 'S': seq[ i ] = 'S'; break;
 114             case 'y': seq[ i ] = 'r'; break;
 115             case 'Y': seq[ i ] = 'R'; break;
 116             case 'k': seq[ i ] = 'm'; break;
 117             case 'K': seq[ i ] = 'M'; break;
 118             case 'b': seq[ i ] = 'v'; break;
 119             case 'B': seq[ i ] = 'V'; break;
 120             case 'd': seq[ i ] = 'h'; break;
 121             case 'D': seq[ i ] = 'H'; break;
 122             case 'h': seq[ i ] = 'd'; break;
 123             case 'H': seq[ i ] = 'D'; break;
 124             case 'v': seq[ i ] = 'b'; break;
 125             case 'V': seq[ i ] = 'B'; break;
 126             case 'n': seq[ i ] = 'n'; break;
 127             case 'N': seq[ i ] = 'N'; break;
 128             default: break;
 129         }
 130     }
 131 }
 132
 133
 134 void complement_rna( char *seq )
 135 {
 136     /* Martin A. Hansen, May 2008 */
 137
 138     /* Complements an RNA sequence including */
 139     /* ambiguity coded nucleotides. */;
 140
 141     size_t i;
 142
 143     for ( i = 0; seq[ i ]; i++ )
 144     {
 145         switch ( seq[ i ] )
 146         {
 147             case 'a': seq[ i ] = 'u'; break;
 148             case 'A': seq[ i ] = 'U'; break;
 149             case 'c': seq[ i ] = 'g'; break;
 150             case 'C': seq[ i ] = 'G'; break;
 151             case 'g': seq[ i ] = 'c'; break;
 152             case 'G': seq[ i ] = 'C'; break;
 153             case 't': seq[ i ] = 'a'; break;
 154             case 'u': seq[ i ] = 'a'; break;
 155             case 'T': seq[ i ] = 'A'; break;
 156             case 'U': seq[ i ] = 'A'; break;
 157             case 'm': seq[ i ] = 'k'; break;
 158             case 'M': seq[ i ] = 'K'; break;
 159             case 'r': seq[ i ] = 'y'; break;
 160             case 'R': seq[ i ] = 'Y'; break;
 161             case 'w': seq[ i ] = 'w'; break;
 162             case 'W': seq[ i ] = 'W'; break;
 163             case 's': seq[ i ] = 'S'; break;
 164             case 'S': seq[ i ] = 'S'; break;
 165             case 'y': seq[ i ] = 'r'; break;
 166             case 'Y': seq[ i ] = 'R'; break;
 167             case 'k': seq[ i ] = 'm'; break;
 168             case 'K': seq[ i ] = 'M'; break;
 169             case 'b': seq[ i ] = 'v'; break;
 170             case 'B': seq[ i ] = 'V'; break;
 171             case 'd': seq[ i ] = 'h'; break;
 172             case 'D': seq[ i ] = 'H'; break;
 173             case 'h': seq[ i ] = 'd'; break;
 174             case 'H': seq[ i ] = 'D'; break;
 175             case 'v': seq[ i ] = 'b'; break;
 176             case 'V': seq[ i ] = 'B'; break;
 177             case 'n': seq[ i ] = 'n'; break;
 178             case 'N': seq[ i ] = 'N'; break;
 179             default: break;
 180         }
 181     }
 182 }
 183
 184
 185 void reverse( char *string )
 186 {
 187     /* Martin A. Hansen, May 2008 */
 188
 189     /* Reverses a string in place. */
 190
 191     char   c;
 192     size_t i;
 193     size_t j;
 194
 195     i = 0;
 196     j = strlen( string ) - 1;
 197
 198     while ( i <= j )
 199     {
 200         c = string[ i ];
 201
 202         string[ i ] = string[ j ];
 203         string[ j ] = c;
 204
 205         i++;
 206         j--;
 207     }
 208 }
 209
 210
 211 void seq2nuc_simple( char *seq )
 212 {
 213     /* Martin A. Hansen, May 2008 */
 214
 215     /* Uppercases all DNA letters, while transforming */
 216     /* all non-DNA letters in sequence to Ns. */
 217
 218     size_t i;
 219
 220     for ( i = 0; seq[ i ]; i++ )
 221     {
 222         switch ( seq[ i ] )
 223         {
 224             case 'A': break;
 225             case 'T': break;
 226             case 'C': break;
 227             case 'G': break;
 228             case 'U': break;
 229             case 'N': break;
 230             case 'a': seq[ i ] = 'A'; break;
 231             case 't': seq[ i ] = 'T'; break;
 232             case 'c': seq[ i ] = 'C'; break;
 233             case 'g': seq[ i ] = 'G'; break;
 234             case 'u': seq[ i ] = 'U'; break;
 235             default:  seq[ i ] = 'N';
 236         }
 237     }
 238 }
 239
 240
 241 void dna2rna( char *seq )
 242 {
 243     /* Martin A. Hansen, May 2008 */
 244
 245     /* Converts a DNA sequence to RNA by changing T and t to U and u. */
 246
 247     size_t i;
 248
 249     for ( i = 0; seq[ i ]; i++ )
 250     {
 251         switch ( seq[ i ] )
 252         {
 253             case 't': seq[ i ] = 'u'; break;
 254             case 'T': seq[ i ] = 'U'; break;
 255             default: break;
 256         }
 257     }
 258 }
 259
 260
 261 void rna2dna( char *seq )
 262 {
 263     /* Martin A. Hansen, May 2008 */
 264
 265     /* Converts a RNA sequence to RNA by changing T and u to T and t. */
 266
 267     size_t i;
 268
 269     for ( i = 0; seq[ i ]; i++ )
 270     {
 271         switch ( seq[ i ] )
 272         {
 273             case 'u': seq[ i ] = 't'; break;
 274             case 'U': seq[ i ] = 'T'; break;
 275             default: break;
 276         }
 277     }
 278 }
 279
 280
 281 bool is_dna( char *seq )
 282 {
 283     /* Martin A. Hansen, May 2008 */
 284
 285     /* Determines if a given sequence is DNA, */
 286     /* from inspection of the first 100 residues. */
 287
 288     size_t i;
 289
 290     for ( i = 0; seq[ i ]; i++ )
 291     {
 292         switch ( seq[ i ] )
 293         {
 294             case 'A': case 'a': break;
 295             case 'G': case 'g': break;
 296             case 'C': case 'c': break;
 297             case 'T': case 't': break;
 298             case 'R': case 'r': break;
 299             case 'Y': case 'y': break;
 300             case 'W': case 'w': break;
 301             case 'S': case 's': break;
 302             case 'M': case 'm': break;
 303             case 'K': case 'k': break;
 304             case 'H': case 'h': break;
 305             case 'D': case 'd': break;
 306             case 'V': case 'v': break;
 307             case 'B': case 'b': break;
 308             case 'N': case 'n': break;
 309             case '-': break;
 310             case '~': break;
 311             case '_': break;
 312             case '.': break;
 313             default: return FALSE;
 314         }
 315
 316         if ( i == 100 ) {
 317             break;
 318         }
 319     }
 320
 321     return TRUE;
 322 }
 323
 324
 325 bool is_rna( char *seq )
 326 {
 327     /* Martin A. Hansen, May 2008 */
 328
 329     /* Determines if a given sequence is RNA, */
 330     /* from inspection of the first 100 residues. */
 331
 332     size_t i;
 333
 334     for ( i = 0; seq[ i ]; i++ )
 335     {
 336         switch ( seq[ i ] )
 337         {
 338             case 'A': case 'a': break;
 339             case 'G': case 'g': break;
 340             case 'C': case 'c': break;
 341             case 'U': case 'u': break;
 342             case 'R': case 'r': break;
 343             case 'Y': case 'y': break;
 344             case 'W': case 'w': break;
 345             case 'S': case 's': break;
 346             case 'M': case 'm': break;
 347             case 'K': case 'k': break;
 348             case 'H': case 'h': break;
 349             case 'D': case 'd': break;
 350             case 'V': case 'v': break;
 351             case 'B': case 'b': break;
 352             case 'N': case 'n': break;
 353             case '-': break;
 354             case '~': break;
 355             case '_': break;
 356             case '.': break;
 357             default: return FALSE;
 358         }
 359
 360         if ( i == 100 ) {
 361             break;
 362         }
 363     }
 364
 365     return TRUE;
 366 }
 367
 368
 369 bool is_protein( char *seq )
 370 {
 371     /* Martin A. Hansen, May 2008 */
 372
 373     /* Determines if a given sequence is protein, */
 374     /* from inspection of the first 100 residues. */
 375
 376     size_t i;
 377
 378     for ( i = 0; seq[ i ]; i++ )
 379     {
 380         switch ( seq[ i ] )
 381         {
 382             case 'K': case 'k': break;
 383             case 'R': case 'r': break;
 384             case 'H': case 'h': break;
 385             case 'D': case 'd': break;
 386             case 'E': case 'e': break;
 387             case 'S': case 's': break;
 388             case 'T': case 't': break;
 389             case 'N': case 'n': break;
 390             case 'Q': case 'q': break;
 391             case 'A': case 'a': break;
 392             case 'V': case 'v': break;
 393             case 'I': case 'i': break;
 394             case 'L': case 'l': break;
 395             case 'M': case 'm': break;
 396             case 'F': case 'f': break;
 397             case 'Y': case 'y': break;
 398             case 'W': case 'w': break;
 399             case 'C': case 'c': break;
 400             case 'G': case 'g': break;
 401             case 'P': case 'p': break;
 402             case 'Z': case 'z': break;
 403             case 'B': case 'b': break;
 404             case 'X': case 'x': break;
 405             case '*': break;
 406             case '-': break;
 407             case '~': break;
 408             case '_': break;
 409             case '.': break;
 410             default: return FALSE;
 411         }
 412
 413         if ( i == 100 ) {
 414             break;
 415         }
 416     }
 417
 418     return TRUE;
 419 }
 420
 421
 422 char *seq_guess_type( char *seq )
 423 {
 424     /* Martin A. Hansen, May 2008 */
 425
 426     /* Guess the type of a given sequnce, */
 427     /* which is returned as a pointer to a string. */
 428
 429     char *type;
 430
 431     type = mem_get( 8 );
 432
 433     if ( is_dna( seq ) ) {
 434         type = "DNA";
 435     } else if ( is_rna( seq ) ) {
 436         type = "RNA";
 437     } else if ( is_protein( seq ) ) {
 438         type = "PROTEIN";
 439     } else {
 440         die( "Could not guess sequence type.\n" );
 441     }
 442
 443     return type;
 444 }
 445
 446
 447 bool contain_N( char *seq )
 448 {
 449     /* Martin A. Hansen, May 2008 */
 450
 451     /* Check if a sequence contain N or n residues. */
 452
 453     size_t i;
 454
 455     for ( i = 0; seq[ i ]; i++ )
 456     {
 457         switch ( seq[ i ] )
 458         {
 459             case 'N': case 'n': return TRUE;
 460             default: break;
 461         }
 462     }
 463
 464     return FALSE;
 465 }
 466
 467
 468 int oligo2bin( char *oligo )
 469 {
 470     /* Martin A. Hansen, August 2004 */
 471
 472     /* Pack a max 15 nucleotide long oligo into a four byte integer. */
 473
 474     int i;
 475     int bin;
 476
 477     if ( strlen( oligo ) > 15 ) {
 478         die( "Oligo will not fit in an integer." );
 479     }
 480
 481     bin = 0;
 482
 483     for ( i = 0; oligo[ i ]; i++ )
 484     {
 485         bin <<= 2;
 486
 487         switch ( oligo[ i ] )
 488         {
 489             case 'A': case 'a': bin |= 0; break;
 490             case 'N': case 'n': bin |= 0; break;
 491             case 'T': case 't': bin |= 1; break;
 492             case 'U': case 'u': bin |= 1; break;
 493             case 'C': case 'c': bin |= 2; break;
 494             case 'G': case 'g': bin |= 3; break;
 495             default: die( "Unrecognized nucleotide." );
 496         }
 497     }
 498
 499     return bin;
 500 }