about summary refs log tree commit diff
path: root/src/mathfunc.cpp
diff options
context:
space:
mode:
Diffstat (limited to 'src/mathfunc.cpp')
-rw-r--r--src/mathfunc.cpp33
1 files changed, 33 insertions, 0 deletions
diff --git a/src/mathfunc.cpp b/src/mathfunc.cpp
index a351509..c82872e 100644
--- a/src/mathfunc.cpp
+++ b/src/mathfunc.cpp
@@ -645,3 +645,36 @@ double UcharToDouble02(const unsigned char c) { return (double)c * 0.01; }
 unsigned char Double02ToUchar(const double dosage) {
   return (int)(dosage * 100);
 }
+
+
+std::tuple<size_t, size_t, size_t> compute_ratio(size_t ni_total, const gsl_vector *gs) {
+  size_t n_0 = 0;
+  size_t n_1 = 0;
+  size_t n_2 = 0;
+  for (size_t i = 0; i < ni_total; ++i) { // read genotypes
+    double geno = gsl_vector_get(gs, i);
+    if (geno >= 0 && geno <= 0.5) {
+      n_0++;
+    }
+    if (geno > 0.5 && geno < 1.5) {
+      n_1++;
+    }
+    if (geno >= 1.5 && geno <= 2.0) {
+      n_2++;
+    }
+  }
+  return {n_0,n_1,n_2};
+}
+
+double compute_maf(size_t ni_total, size_t ni_test, size_t n_miss, const double *gs, const vector<int> &indicator) {
+  double maf = 0.0;
+
+  for (size_t i = 0; i < ni_total; ++i) { // read genotypes
+    if (indicator[i]) {
+      double geno = gs[i];
+      maf += geno; // 0..2 range
+    }
+  }
+  maf /= 2.0 * (double)(ni_test - n_miss); // Assumption is that geno value is between 0 and 2. FIXME
+  return maf;
+}