<commit_before>#include "Variant.h"
#include "split.h"
#include "convert.h"
#include <getopt.h>

using namespace std;
using namespace vcf;

bool isTransition(const string& ref, const string& alt) {
    if (((ref == "A" && alt == "G") || (ref == "G" && alt == "A")) ||
        ((ref == "C" && alt == "T") || (ref == "T" && alt == "C"))) {
        return true;
    } else {
        return false;
    }
}

bool isDeamination(const string& ref, const string& alt) {
    if ((ref == "G" && alt == "A") ||
        (ref == "C" && alt == "T")) {
        return true;
    } else {
        return false;
    }
}

bool isAmination(const string& ref, const string& alt) {
    if ((ref == "A" && alt == "G") ||
        (ref == "T" && alt == "C")) {
        return true;
    } else {
        return false;
    }
}

void printSummary(char** argv) {
    cerr << "usage: " << argv[0] << " [options] <vcf file>" << endl
         << endl
         << "    -r, --region          specify a region on which to target the stats, requires a BGZF" << endl
         << "                          compressed file which has been indexed with tabix.  any number of" << endl
         << "                          regions may be specified." << endl
         << endl
         << "Prints statistics about variants in the input VCF file." << endl;
}


int main(int argc, char** argv) {

    vector<string> regions;

    int c;
    while (true) {
        static struct option long_options[] =
        {
            /* These options set a flag. */
            //{"verbose", no_argument,       &verbose_flag, 1},
            {"help", no_argument, 0, 'h'},
            {"region", required_argument, 0, 'r'},
            //{"length",  no_argument, &printLength, true},
            {0, 0, 0, 0}
        };
        /* getopt_long stores the option index here. */
        int option_index = 0;

        c = getopt_long (argc, argv, "hr:",
                         long_options, &option_index);

      /* Detect the end of the options. */
          if (c == -1)
            break;
 
          switch (c)
            {
            case 0:
            /* If this option set a flag, do nothing else now. */
            if (long_options[option_index].flag != 0)
              break;
            printf ("option %s", long_options[option_index].name);
            if (optarg)
              printf (" with arg %s", optarg);
            printf ("\n");
            break;

          case 'h':
            printSummary(argv);
            exit(0);
            break;

          case 'r':
            regions.push_back(optarg);
            break;
          
          default:
            abort ();
          }
      }

    VariantCallFile variantFile;
    string inputFilename;
    if (optind == argc - 1) {
        inputFilename = argv[optind];
        variantFile.open(inputFilename);
    } else {
        variantFile.open(std::cin);
    }

    if (!variantFile.is_open()) {
        return 1;
    }

    Variant var(variantFile);

    vector<string>::iterator regionItr = regions.begin();

    int variantAlleles = 0;
    int uniqueVariantAlleles = 0;
    int variantSites = 0;
    int snps = 0;
    int transitions = 0;
    int transversions = 0;
    int deaminations = 0;
    int aminations = 0;
    int totalinsertions = 0;
    int totaldeletions = 0;
    int insertedbases = 0;
    int deletedbases = 0;
    int totalmnps = 0;
    int totalcomplex = 0;
    int mismatchbases = 0;
    int mnpbases = 0;
    map<int, int> insertions;
    map<int, int> deletions;
    map<int, int> mnps;
    map<int, int> complexsubs;

    do {

        if (!inputFilename.empty() && !regions.empty()) {
            string regionStr = *regionItr++;
            variantFile.setRegion(regionStr);
        }

        while (variantFile.getNextVariant(var)) {
            ++variantSites;
            map<string, vector<VariantAllele> > alternates = var.parsedAlternates();
            map<VariantAllele, vector<string> > uniqueVariants;
            for (vector<string>::iterator a = var.alt.begin(); a != var.alt.end(); ++a) {
                string& alternate = *a;
                vector<VariantAllele>& vav = alternates[alternate];
                if (vav.size() > 1) {
                    ++totalcomplex;
                }
                for (vector<VariantAllele>::iterator v = vav.begin(); v != vav.end(); ++v) {
                    uniqueVariants[*v].push_back(alternate);
                }
            }
            uniqueVariantAlleles += uniqueVariants.size();
            variantAlleles += var.alt.size();

            for (map<VariantAllele, vector<string> >::iterator v = uniqueVariants.begin(); v != uniqueVariants.end(); ++v) {
                const VariantAllele& va = v->first;
                //cout << va.ref << "/" << va.alt << endl;

                if (va.ref.size() == va.alt.size()) {
                    if (va.ref.size() == 1) {
                        ++snps;
                        ++mismatchbases;
                        if (isTransition(va.ref, va.alt)) {
                            ++transitions;
                        } else {
                            ++transversions;
                        }
                        if (isAmination(va.ref, va.alt)) {
                            ++aminations;
                        }
                        if (isDeamination(va.ref, va.alt)) {
                            ++deaminations;
                        }
                    } else {
                        ++totalmnps;
                        ++mnps[va.alt.size()]; // not entirely correct
                        string::const_iterator r = va.ref.begin();
                        for (string::const_iterator a = va.alt.begin(); a != va.alt.end(); ++a, ++r) {
                            string rstr = string(*r, 1);
                            string astr = string(*a, 1);
                            if (rstr == astr) {
                                continue;
                            }
                            if (isTransition(rstr, astr)) {
                                ++transitions;
                            } else {
                                ++transversions;
                            }
                            if (isAmination(rstr, astr)) {
                                ++aminations;
                            }
                            if (isDeamination(rstr, astr)) {
                                ++deaminations;
                            }
                            ++mismatchbases;
                            ++mnpbases;
                        }
                    }
                } else if (va.ref.size() > va.alt.size()) {
                    int diff = va.ref.size() - va.alt.size();
                    deletedbases += diff;
                    ++totaldeletions;
                    ++deletions[diff];
                } else {
                    int diff = va.alt.size() - va.ref.size();
                    insertedbases += diff;
                    ++totalinsertions;
                    ++insertions[diff];
                }<commit_msg>mnp ts/tv and deamination bug<commit_after>#include "Variant.h"
#include "split.h"
#include "convert.h"
#include <getopt.h>

using namespace std;
using namespace vcf;

bool isTransition(const string& ref, const string& alt) {
    if (((ref == "A" && alt == "G") || (ref == "G" && alt == "A")) ||
        ((ref == "C" && alt == "T") || (ref == "T" && alt == "C"))) {
        return true;
    } else {
        return false;
    }
}

bool isDeamination(const string& ref, const string& alt) {
    if ((ref == "G" && alt == "A") ||
        (ref == "C" && alt == "T")) {
        return true;
    } else {
        return false;
    }
}

bool isAmination(const string& ref, const string& alt) {
    if ((ref == "A" && alt == "G") ||
        (ref == "T" && alt == "C")) {
        return true;
    } else {
        return false;
    }
}

void printSummary(char** argv) {
    cerr << "usage: " << argv[0] << " [options] <vcf file>" << endl
         << endl
         << "    -r, --region          specify a region on which to target the stats, requires a BGZF" << endl
         << "                          compressed file which has been indexed with tabix.  any number of" << endl
         << "                          regions may be specified." << endl
         << endl
         << "Prints statistics about variants in the input VCF file." << endl;
}


int main(int argc, char** argv) {

    vector<string> regions;

    int c;
    while (true) {
        static struct option long_options[] =
        {
            /* These options set a flag. */
            //{"verbose", no_argument,       &verbose_flag, 1},
            {"help", no_argument, 0, 'h'},
            {"region", required_argument, 0, 'r'},
            //{"length",  no_argument, &printLength, true},
            {0, 0, 0, 0}
        };
        /* getopt_long stores the option index here. */
        int option_index = 0;

        c = getopt_long (argc, argv, "hr:",
                         long_options, &option_index);

      /* Detect the end of the options. */
          if (c == -1)
            break;
 
          switch (c)
            {
            case 0:
            /* If this option set a flag, do nothing else now. */
            if (long_options[option_index].flag != 0)
              break;
            printf ("option %s", long_options[option_index].name);
            if (optarg)
              printf (" with arg %s", optarg);
            printf ("\n");
            break;

          case 'h':
            printSummary(argv);
            exit(0);
            break;

          case 'r':
            regions.push_back(optarg);
            break;
          
          default:
            abort ();
          }
      }

    VariantCallFile variantFile;
    string inputFilename;
    if (optind == argc - 1) {
        inputFilename = argv[optind];
        variantFile.open(inputFilename);
    } else {
        variantFile.open(std::cin);
    }

    if (!variantFile.is_open()) {
        return 1;
    }

    Variant var(variantFile);

    vector<string>::iterator regionItr = regions.begin();

    int variantAlleles = 0;
    int uniqueVariantAlleles = 0;
    int variantSites = 0;
    int snps = 0;
    int transitions = 0;
    int transversions = 0;
    int deaminations = 0;
    int aminations = 0;
    int totalinsertions = 0;
    int totaldeletions = 0;
    int insertedbases = 0;
    int deletedbases = 0;
    int totalmnps = 0;
    int totalcomplex = 0;
    int mismatchbases = 0;
    int mnpbases = 0;
    map<int, int> insertions;
    map<int, int> deletions;
    map<int, int> mnps;
    map<int, int> complexsubs;

    do {

        if (!inputFilename.empty() && !regions.empty()) {
            string regionStr = *regionItr++;
            variantFile.setRegion(regionStr);
        }

        while (variantFile.getNextVariant(var)) {
            ++variantSites;
            map<string, vector<VariantAllele> > alternates = var.parsedAlternates();
            map<VariantAllele, vector<string> > uniqueVariants;
            for (vector<string>::iterator a = var.alt.begin(); a != var.alt.end(); ++a) {
                string& alternate = *a;
                vector<VariantAllele>& vav = alternates[alternate];
                if (vav.size() > 1) {
                    ++totalcomplex;
                }
                for (vector<VariantAllele>::iterator v = vav.begin(); v != vav.end(); ++v) {
                    uniqueVariants[*v].push_back(alternate);
                }
            }
            uniqueVariantAlleles += uniqueVariants.size();
            variantAlleles += var.alt.size();

            for (map<VariantAllele, vector<string> >::iterator v = uniqueVariants.begin(); v != uniqueVariants.end(); ++v) {
                const VariantAllele& va = v->first;
                //cout << va.ref << "/" << va.alt << endl;

                if (va.ref.size() == va.alt.size()) {
                    if (va.ref.size() == 1) {
                        ++snps;
                        ++mismatchbases;
                        if (isTransition(va.ref, va.alt)) {
                            ++transitions;
                        } else {
                            ++transversions;
                        }
                        if (isAmination(va.ref, va.alt)) {
                            ++aminations;
                        }
                        if (isDeamination(va.ref, va.alt)) {
                            ++deaminations;
                        }
                    } else {
                        ++totalmnps;
                        ++mnps[va.alt.size()]; // not entirely correct
                        string::const_iterator r = va.ref.begin();
                        for (string::const_iterator a = va.alt.begin(); a != va.alt.end(); ++a, ++r) {
                            string rstr = string(1, *r);
                            string astr = string(1, *a);
                            if (rstr == astr) {
                                continue;
                            }
                            if (isTransition(rstr, astr)) {
                                ++transitions;
                            } else {
                                ++transversions;
                            }
                            if (isAmination(rstr, astr)) {
                                ++aminations;
                            }
                            if (isDeamination(rstr, astr)) {
                                ++deaminations;
                            }
                            ++mismatchbases;
                            ++mnpbases;
                        }
                    }
                } else if (va.ref.size() > va.alt.size()) {
                    int diff = va.ref.size() - va.alt.size();
                    deletedbases += diff;
                    ++totaldeletions;
                    ++deletions[diff];
                } else {
                    int diff = va.alt.size() - va.ref.size();
                    insertedbases += diff;
                    ++totalinsertions;
                    ++insertions[diff];
                }
{'fraction_non_alphanumeric': 0.12971331732806102, 'fraction_numerical': 0.004942804688603305, 'mean_word_length': 1.2102059925093633, 'pattern_counts': {'":': 0, '<': 109, '<?xml version=': 0, '>': 43, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 4, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'dir': 'git-commits-cleaned', 'id': '6089254', 'n_tokens_mistral': 4027, 'n_tokens_neox': 3823, 'n_words': 968}