testSegment.cpp 15.5 KB
Newer Older
Mathieu Giraud's avatar
Mathieu Giraud committed
1 2 3 4

#include <fstream>
#include <iostream>
#include <string>
5
#include "core/germline.h"
Mathieu Giraud's avatar
Mathieu Giraud committed
6 7 8 9 10 11 12 13
#include "core/kmerstore.h"
#include "core/dynprog.h"
#include "core/fasta.h"
#include "core/segment.h"
#include "core/windowExtractor.h"

using namespace std;

Mikael Salson's avatar
Mikael Salson committed
14
void testFineSegment()
Mathieu Giraud's avatar
Mathieu Giraud committed
15 16 17 18 19 20 21
{
  Fasta seqV("../../germline/IGHV.fa", 2);
  Fasta seqD("../../germline/IGHD.fa", 2);
  Fasta seqJ("../../germline/IGHJ.fa", 2);
  
  Fasta data("../../data/Stanford_S22.fasta", 1, " ");

22
  Germline *germline ;
23 24
  germline = new Germline("IGH", 'G', seqV, seqD, seqJ, 0, 50);
  germline->new_index("####");
25

Mathieu Giraud's avatar
Mathieu Giraud committed
26 27 28
  Sequence seq = data.read(2);
      
  //segmentation VJ
29
  FineSegmenter s(seq, germline, VDJ);
Mathieu Giraud's avatar
Mathieu Giraud committed
30 31 32 33
	
  TAP_TEST(s.isSegmented(), TEST_SEGMENT_POSITION, "is segmented (VJ)") ;
  
  //segmentation D
34
  s.FineSegmentD(germline);
Mathieu Giraud's avatar
Mathieu Giraud committed
35 36 37 38 39
  
  TAP_TEST(s.isSegmented(), TEST_SEGMENT_POSITION, "is segmented (VDJ)") ;

  // Revcomp sequence and tests that the results are the same.
  seq.sequence = revcomp(seq.sequence);
40
  FineSegmenter s2(seq, germline, VDJ);
Mathieu Giraud's avatar
Mathieu Giraud committed
41 42 43

  TAP_TEST(s2.isSegmented(), TEST_SEGMENT_POSITION, "is segmented (VJ)") ;
  //segmentation D
44
  s2.FineSegmentD(germline);
Mathieu Giraud's avatar
Mathieu Giraud committed
45 46 47 48 49 50 51 52 53 54 55 56 57 58 59
  TAP_TEST(s2.isSegmented(), TEST_SEGMENT_POSITION, "is segmented (VDJ)") ;
  
  TAP_TEST(s.getLeft() == s2.getLeft(), TEST_SEGMENT_REVCOMP, " left segmentation position");
  TAP_TEST(s.getRight() == s2.getRight(), TEST_SEGMENT_REVCOMP, " right segmentation position");
  TAP_TEST(s.getLeftD() == s2.getLeftD(), TEST_SEGMENT_REVCOMP, " left D segmentation position");
  TAP_TEST(s.getRightD() == s2.getRightD(), TEST_SEGMENT_REVCOMP, " right D segmentation position");
  TAP_TEST(s.isReverse() == !s2.isReverse(), TEST_SEGMENT_REVCOMP, " sequence reversed");
  TAP_TEST(s.info.substr(1) == s2.info.substr(1), TEST_SEGMENT_REVCOMP, " info string " << endl <<
           "s  = " << s.info << endl <<
	   "s2 = " << s2.info);
  TAP_TEST(s.info.substr(0,1) != s2.info.substr(0,1), TEST_SEGMENT_REVCOMP, "first character (strand) of info string " << endl <<
           "s  = " << s.info << endl <<
	   "s2 = " << s2.info);

  TAP_TEST(s.code == s2.code, TEST_SEGMENT_REVCOMP, " code string");
60
  delete germline;
Mathieu Giraud's avatar
Mathieu Giraud committed
61 62 63 64 65 66 67 68 69 70 71 72
}

/**
 * Test segmentation when there is an overlap between V and J (and no N)
 */
void testSegmentOverlap()
{
  Fasta seqV("../../germline/TRGV.fa", 2);
  Fasta seqJ("../../germline/TRGJ.fa", 2);
  
  Fasta data("../../data/bug-segment-overlap.fa", 1, " ");
  
73
  Germline *germline1 ;
74 75 76
  germline1 = new Germline("TRG", 'G', seqV, seqV, seqJ, -50, 50);
  germline1->new_index("##########");

77
  Germline *germline2 ;
78 79
  germline2 = new Germline("TRG2", 'G', seqV, seqV, seqJ, -50, 50);
  germline2->new_index("##########");
80

Mathieu Giraud's avatar
Mathieu Giraud committed
81
  for (int i = 0; i < data.size(); i++) {
82 83
    KmerSegmenter ks(data.read(i), germline1);

Mathieu Giraud's avatar
Mathieu Giraud committed
84 85 86 87 88
    TAP_TEST(ks.seg_V + ks.seg_N + ks.seg_J == data.sequence(i)
             || ks.seg_V + ks.seg_N + ks.seg_J == revcomp(data.sequence(i)), 
             TEST_KMER_SEGMENT_OVERLAP,
             " V= " << ks.seg_V << ", N = " << ks.seg_N << ", J = " << ks.seg_J);

89
    FineSegmenter fs(data.read(i), germline2, VDJ); 
Mathieu Giraud's avatar
Mathieu Giraud committed
90 91 92 93 94
    TAP_TEST(fs.seg_V + fs.seg_N + fs.seg_J == data.sequence(i)
             || fs.seg_V + fs.seg_N + fs.seg_J == revcomp(data.sequence(i)), 
             TEST_FINE_SEGMENT_OVERLAP,
             " V= " << fs.seg_V << ", N = " << fs.seg_N << ", J = " << fs.seg_J);
  }
95

96
  delete germline1;
97
  delete germline2;
Mathieu Giraud's avatar
Mathieu Giraud committed
98 99 100 101 102 103 104 105
}

void testSegmentationCause() {
  Fasta seqV("../../germline/TRGV.fa", 2);
  Fasta seqJ("../../germline/TRGJ.fa", 2);
  
  Fasta data("../../data/segmentation.fasta", 1, " ");

106
  Germline *germline ;
107 108
  germline = new Germline("TRG", 'G', seqV, seqV, seqJ, 0, 10);
  germline->new_index("##########");
109

Mathieu Giraud's avatar
Mathieu Giraud committed
110 111 112
  int nb_checked = 0;

  for (int i = 0; i < data.size(); i++) {
113 114
    KmerSegmenter ks(data.read(i), germline);
    
Mathieu Giraud's avatar
Mathieu Giraud committed
115
    if (data.read(i).label == "seq-seg+") {
116
      TAP_TEST(ks.isSegmented(), TEST_KMER_IS_SEGMENTED, "seq is " << data.label(i));
117
      TAP_TEST(ks.getSegmentationStatus() == SEG_PLUS, TEST_KMER_SEGMENTATION_CAUSE, ks.getInfoLineWithAffects());
118 119
      TAP_TEST(ks.getJunction(30) == "GCCACCTGGGACAGGGAATTATTATAAGAA"
               || ks.getJunction(30) == "TGCCACCTGGGACAGGGAATTATTATAAGA", 
120
               TEST_KMER_JUNCTION, ks.getInfoLineWithAffects());
121 122
      TAP_TEST(ks.getLeft() == 17, TEST_KMER_LEFT, "left = " << ks.getLeft());
      TAP_TEST(ks.getRight() == 18, TEST_KMER_RIGHT, "right = " << ks.getRight());
123

124
      ks.setSegmentationStatus(NOT_PROCESSED);
125 126
      TAP_TEST(! ks.isSegmented(), TEST_SET_SEGMENTATION_CAUSE, ks.getInfoLineWithAffects());
      TAP_TEST(ks.getSegmentationStatus() == NOT_PROCESSED, TEST_SET_SEGMENTATION_CAUSE, ks.getInfoLineWithAffects());
127
      ks.setSegmentationStatus(UNSEG_NOISY);
128 129
      TAP_TEST(! ks.isSegmented(), TEST_SET_SEGMENTATION_CAUSE, ks.getInfoLineWithAffects());
      TAP_TEST(ks.getSegmentationStatus() == UNSEG_NOISY, TEST_SET_SEGMENTATION_CAUSE, ks.getInfoLineWithAffects());
130
      ks.setSegmentationStatus(SEG_PLUS);
131 132
      TAP_TEST(ks.isSegmented(), TEST_SET_SEGMENTATION_CAUSE, ks.getInfoLineWithAffects());
      TAP_TEST(ks.getSegmentationStatus(), TEST_SET_SEGMENTATION_CAUSE, ks.getInfoLineWithAffects());
Mathieu Giraud's avatar
Mathieu Giraud committed
133 134
      nb_checked++;
    } else if (data.read(i).label == "seq-seg-") {
135 136
      TAP_TEST(ks.isSegmented(), TEST_KMER_IS_SEGMENTED, ks.getInfoLineWithAffects());
      TAP_TEST(ks.getSegmentationStatus() == SEG_MINUS, TEST_KMER_SEGMENTATION_CAUSE, ks.getInfoLineWithAffects());
137 138
      TAP_TEST(ks.getJunction(30) == "GCCACCTGGGACAGGGAATTATTATAAGAA"
               || ks.getJunction(30) == "TGCCACCTGGGACAGGGAATTATTATAAGA", 
139
               TEST_KMER_JUNCTION, ks.getInfoLineWithAffects());
140 141
      TAP_TEST(ks.getLeft() == 17, TEST_KMER_LEFT, "left = " << ks.getLeft());
      TAP_TEST(ks.getRight() == 18, TEST_KMER_RIGHT, "right = " << ks.getRight());
Mathieu Giraud's avatar
Mathieu Giraud committed
142 143
      nb_checked++;
    } else if (data.read(i).label == "seq-short") {
144 145
      TAP_TEST(! ks.isSegmented(), TEST_KMER_IS_SEGMENTED, ks.getInfoLineWithAffects());
      TAP_TEST(ks.getSegmentationStatus() == UNSEG_TOO_SHORT, TEST_KMER_SEGMENTATION_CAUSE, ks.getInfoLineWithAffects());
Mathieu Giraud's avatar
Mathieu Giraud committed
146 147
      nb_checked++;
    } else if (data.read(i).label == "seq-strand") {
148 149
      TAP_TEST(! ks.isSegmented(), TEST_KMER_IS_SEGMENTED, ks.getInfoLineWithAffects());
      TAP_TEST(ks.getSegmentationStatus() == UNSEG_STRAND_NOT_CONSISTENT, TEST_KMER_SEGMENTATION_CAUSE, ks.getInfoLineWithAffects());
Mathieu Giraud's avatar
Mathieu Giraud committed
150 151
      nb_checked++;
    } else if (data.read(i).label == "seq-zero") {
152 153
      TAP_TEST(! ks.isSegmented(), TEST_KMER_IS_SEGMENTED, ks.getInfoLineWithAffects());
      TAP_TEST(ks.getSegmentationStatus() == UNSEG_TOO_FEW_ZERO, TEST_KMER_SEGMENTATION_CAUSE, ks.getInfoLineWithAffects());
Mathieu Giraud's avatar
Mathieu Giraud committed
154 155
      nb_checked++;
    } else if (data.read(i).label == "seq-fewV") {
156 157
      TAP_TEST(! ks.isSegmented(), TEST_KMER_IS_SEGMENTED, ks.getInfoLineWithAffects());
      TAP_TEST(ks.getSegmentationStatus() == UNSEG_TOO_FEW_V, TEST_KMER_SEGMENTATION_CAUSE, ks.getInfoLineWithAffects());
Mathieu Giraud's avatar
Mathieu Giraud committed
158 159
      nb_checked++;
    } else if (data.read(i).label == "seq-fewJ") {
160 161
      TAP_TEST(! ks.isSegmented(), TEST_KMER_IS_SEGMENTED, "fewJ: " << ks.getInfoLineWithAffects());
      TAP_TEST(ks.getSegmentationStatus() == UNSEG_TOO_FEW_J, TEST_KMER_SEGMENTATION_CAUSE, ks.getInfoLineWithAffects());
Mathieu Giraud's avatar
Mathieu Giraud committed
162
      nb_checked++;
163
    } else if (data.read(i).label == "seq-fewV2") {
164 165
      TAP_TEST(! ks.isSegmented(), TEST_KMER_IS_SEGMENTED, ks.getInfoLineWithAffects());
      TAP_TEST(ks.getSegmentationStatus() == UNSEG_TOO_FEW_V, TEST_KMER_SEGMENTATION_CAUSE, ks.getInfoLineWithAffects());
166 167
      nb_checked++;
    } else if (data.read(i).label == "seq-fewJ2") {
168 169
      TAP_TEST(! ks.isSegmented(), TEST_KMER_IS_SEGMENTED, "fewJ: " << ks.getInfoLineWithAffects());
      TAP_TEST(ks.getSegmentationStatus() == UNSEG_TOO_FEW_J, TEST_KMER_SEGMENTATION_CAUSE, ks.getInfoLineWithAffects());
170
      nb_checked++;
171 172 173
    } else if (data.read(i).label == "seq-delta-min-old") {
      // This test was a test for delta_min but with the CountKmerAffectAnalyser
      // the read is segmented, now. So we keep it, but change the test
174 175
      TAP_TEST(ks.isSegmented(), TEST_KMER_IS_SEGMENTED, "delta-min: " << ks.getInfoLineWithAffects());
      TAP_TEST(ks.getSegmentationStatus() == SEG_PLUS, TEST_KMER_SEGMENTATION_CAUSE, ks.getInfoLineWithAffects());
176 177 178
      TAP_TEST(ks.getJunction(30) == "GCCACCTGGGACAGGGAATTATTATAAGAA"
               || ks.getJunction(30) == "TGCCACCTGGGACAGGGAATTATTATAAGA", 
               TEST_KMER_JUNCTION, "junction: " << ks.getJunction(30));
179
      nb_checked++;
Mathieu Giraud's avatar
Mathieu Giraud committed
180
    } else if (data.read(i).label == "seq-delta-min") {
181 182
      TAP_TEST(ks.isSegmented(), TEST_KMER_IS_SEGMENTED, "delta-min: " << ks.getInfoLineWithAffects());
      TAP_TEST(ks.getSegmentationStatus() == SEG_PLUS, TEST_KMER_SEGMENTATION_CAUSE, ks.getInfoLineWithAffects());
183 184
      TAP_TEST(ks.getJunction(21) == "GGCAGTTGGAACAACACTTGT",
               TEST_KMER_JUNCTION, "window: " << ks.getJunction(21));
185 186
      TAP_TEST(ks.getLeft() == 9, TEST_KMER_LEFT, "left = " << ks.getLeft() << ", " << ks.getInfoLineWithAffects());
      TAP_TEST(ks.getRight() == 19, TEST_KMER_RIGHT, "right = " << ks.getRight() << ", " << ks.getInfoLineWithAffects());
Mathieu Giraud's avatar
Mathieu Giraud committed
187 188
      nb_checked++;
    } else if (data.read(i).label == "seq-delta-max") {
189 190
      TAP_TEST(! ks.isSegmented(), TEST_KMER_IS_SEGMENTED, ks.getInfoLineWithAffects());
      TAP_TEST(ks.getSegmentationStatus() == UNSEG_BAD_DELTA_MAX, TEST_KMER_SEGMENTATION_CAUSE, ks.getInfoLineWithAffects());
Mathieu Giraud's avatar
Mathieu Giraud committed
191 192
      nb_checked++;
    } else if (data.read(i).label == "seq-seg-no-window") {
193 194
      TAP_TEST(ks.isSegmented(), TEST_KMER_IS_SEGMENTED, ks.getInfoLineWithAffects());
      TAP_TEST(ks.getSegmentationStatus() == SEG_PLUS, TEST_KMER_SEGMENTATION_CAUSE, ks.getInfoLineWithAffects());
195 196
      TAP_TEST(ks.getLeft() == 11, TEST_KMER_LEFT, "left = " << ks.getLeft());
      TAP_TEST(ks.getRight() == 12, TEST_KMER_RIGHT, "right = " << ks.getRight());
197
      TAP_TEST(ks.getJunction(30) == "", TEST_KMER_JUNCTION, ks.getInfoLineWithAffects());
198 199
      TAP_TEST(ks.getJunction(20) == "CTGGGACAGGGAATTATTAT"
               || ks.getJunction(20) == "CCTGGGACAGGGAATTATTA", TEST_KMER_JUNCTION,"window: " << ks.getJunction(20));
Mathieu Giraud's avatar
Mathieu Giraud committed
200 201 202 203
      nb_checked++;
    }
  }
  
204
  TAP_TEST(nb_checked == 13, TEST_KMER_DATA, "");
205

206
  delete germline;
Mathieu Giraud's avatar
Mathieu Giraud committed
207 208 209 210 211 212 213 214
}

void testExtractor() {
  Fasta seqV("../../germline/TRGV.fa", 2);
  Fasta seqJ("../../germline/TRGJ.fa", 2);
  
  OnlineFasta data("../../data/segmentation.fasta", 1, " ");

215
  Germline *germline ;
216
  germline = new Germline("TRG", 'G', seqV, seqV, seqJ, 0, 10, 0);
217 218
  germline->new_index("##########");

219
  MultiGermline *multi ;
220 221
  multi = new MultiGermline();
  multi->insert(germline);
Mathieu Giraud's avatar
Mathieu Giraud committed
222

223
  WindowExtractor we(multi);
Mathieu Giraud's avatar
Mathieu Giraud committed
224 225 226 227 228 229
  map<string, string> labels;
  ofstream out_seg("segmented.log");
  ofstream out_unseg("unsegmented.log");
  we.setSegmentedOutput(&out_seg);
  we.setUnsegmentedOutput(&out_unseg);

230
  WindowsStorage *ws = we.extract(&data, 30, labels);
231
  // we.out_stats(cout);
Mathieu Giraud's avatar
Mathieu Giraud committed
232

233
  TAP_TEST(we.getNbReads() == 13, TEST_EXTRACTOR_NB_READS, "");
Mathieu Giraud's avatar
Mathieu Giraud committed
234

235
  TAP_TEST(we.getNbSegmented(SEG_PLUS) == 2, TEST_EXTRACTOR_NB_SEGMENTED, "segPlus: " << we.getNbSegmented(SEG_PLUS));
Mathieu Giraud's avatar
Mathieu Giraud committed
236 237 238 239
  TAP_TEST(we.getNbSegmented(SEG_MINUS) == 1, TEST_EXTRACTOR_NB_SEGMENTED, "");
  TAP_TEST(we.getNbSegmented(UNSEG_TOO_SHORT) == 1, TEST_EXTRACTOR_NB_SEGMENTED, "");
  TAP_TEST(we.getNbSegmented(UNSEG_STRAND_NOT_CONSISTENT) == 1, TEST_EXTRACTOR_NB_SEGMENTED, "");
  TAP_TEST(we.getNbSegmented(UNSEG_TOO_FEW_ZERO) == 1, TEST_EXTRACTOR_NB_SEGMENTED, "");
240 241
  TAP_TEST(we.getNbSegmented(UNSEG_TOO_FEW_V) == 2, TEST_EXTRACTOR_NB_SEGMENTED, "");
  TAP_TEST(we.getNbSegmented(UNSEG_TOO_FEW_J) == 2, TEST_EXTRACTOR_NB_SEGMENTED, "");
242
  TAP_TEST(we.getNbSegmented(UNSEG_BAD_DELTA_MIN) == 0, TEST_EXTRACTOR_NB_SEGMENTED, "");
Mathieu Giraud's avatar
Mathieu Giraud committed
243
  TAP_TEST(we.getNbSegmented(UNSEG_BAD_DELTA_MAX) == 1, TEST_EXTRACTOR_NB_SEGMENTED, "");
244
  TAP_TEST(we.getNbSegmented(UNSEG_TOO_SHORT_FOR_WINDOW) == 2, TEST_EXTRACTOR_NB_SEGMENTED, "");
245
  TAP_TEST(we.getNbSegmented(TOTAL_SEG_AND_WINDOW) == 3, TEST_EXTRACTOR_NB_SEGMENTED, "");
Mathieu Giraud's avatar
Mathieu Giraud committed
246

247
  TAP_TEST(we.getAverageSegmentationLength(SEG_PLUS) == 54, TEST_EXTRACTOR_AVG_LENGTH, "average: " << we.getAverageSegmentationLength(SEG_PLUS));
Mathieu Giraud's avatar
Mathieu Giraud committed
248 249 250 251
  TAP_TEST(we.getAverageSegmentationLength(SEG_MINUS) == 36, TEST_EXTRACTOR_AVG_LENGTH, "");
  TAP_TEST(we.getAverageSegmentationLength(UNSEG_TOO_SHORT) == 4, TEST_EXTRACTOR_AVG_LENGTH, "");
  TAP_TEST(we.getAverageSegmentationLength(UNSEG_STRAND_NOT_CONSISTENT) == 36, TEST_EXTRACTOR_AVG_LENGTH, "");
  TAP_TEST(we.getAverageSegmentationLength(UNSEG_TOO_FEW_ZERO) == 36, TEST_EXTRACTOR_AVG_LENGTH, "");
252 253
  TAP_TEST(we.getAverageSegmentationLength(UNSEG_TOO_FEW_V) == 51, TEST_EXTRACTOR_AVG_LENGTH, "");
  TAP_TEST(we.getAverageSegmentationLength(UNSEG_TOO_FEW_J) == 55, TEST_EXTRACTOR_AVG_LENGTH, "");
Mathieu Giraud's avatar
Mathieu Giraud committed
254
  TAP_TEST(we.getAverageSegmentationLength(UNSEG_BAD_DELTA_MAX) == 66, TEST_EXTRACTOR_AVG_LENGTH, "");
255
  TAP_TEST(we.getAverageSegmentationLength(UNSEG_TOO_SHORT_FOR_WINDOW) == 28.5, TEST_EXTRACTOR_AVG_LENGTH, "");
256
  TAP_TEST(we.getAverageSegmentationLength(TOTAL_SEG_AND_WINDOW) == 48, TEST_EXTRACTOR_AVG_LENGTH, "");
Mathieu Giraud's avatar
Mathieu Giraud committed
257 258 259 260 261

  TAP_TEST(out_seg.tellp() > 0, TEST_EXTRACTOR_OUT_SEG, "");
  TAP_TEST(out_unseg.tellp() > 0, TEST_EXTRACTOR_OUT_UNSEG, "");

  delete ws;
262
  delete multi;
Mathieu Giraud's avatar
Mathieu Giraud committed
263
}
Mikael Salson's avatar
Mikael Salson committed
264

265 266 267 268 269 270 271 272 273 274 275 276 277 278 279 280 281 282 283 284 285 286 287 288 289 290
void testProbability() {
  string v_seq[] = {"AAAA", "AAAC", "AAAG", "AAAT", "AACA", "AACC",
                "AACG", "AACT", "AAGA", "AAGC", "AAGG", "AAGT",
                    "AATA", "AATC", "AATG", "AATT", "ACAA", "ACAC",
                "ACAG", "ACAT", "ACCA", "ACCC", "ACCG", "ACCT",
                "ACGA", "ACGC", "ACGG", "ACGT", "ACTA", "ACTC",
                "ACTG", "ACTT", "AGAA", "AGAC", "AGAG", "AGAT",
                "AGCA", "AGCC", "AGCG", "AGCT", "AGGA", "AGGC",
                "AGGG", "AGGT", "AGTA", "AGTC", "AGTG", "AGTT",
                "ATAA", "ATAC", "ATAG", "ATAT", "ATCA", "ATCC",
                "ATCG", "ATCT", "ATGA", "ATGC", "ATGG", "ATGT",
                "ATTA", "ATTC", "ATTG", "ATTT"};
  string j_seq[] = {"CAAA", "CAAC",
                "CAAG", "CAAT", "CACA", "CACC", "CACG", "CACT",
                "CAGA", "CAGC", "CAGG", "CAGT", "CATA", "CATC",
                "CATG", "CATT", "CCAA", "CCAC", "CCAG", "CCAT",
                "CCCA", "CCCC", "CCCG", "CCCT", "CCGA", "CCGC",
                "CCGG", "CCGT", "CCTA", "CCTC", "CCTG", "CCTT",
                "CGAA", "CGAC", "CGAG", "CGAT", "CGCA", "CGCC",
                "CGCG", "CGCT", "CGGA", "CGGC", "CGGG", "CGGT",
                "CGTA", "CGTC", "CGTG", "CGTT", "CTAA", "CTAC",
                "CTAG", "CTAT", "CTCA", "CTCC", "CTCG", "CTCT",
                "CTGA", "CTGC", "CTGG", "CTGT", "CTTA", "CTTC",
                "CTTG", "CTTT"};
  Fasta V, J;
  for (int i = 0; i < 64; i++) {
291 292
    Sequence v = {"V_" + string_of_int(i+33), "V" + string_of_int(i+33), v_seq[i], "", NULL};
    Sequence j = {"J_" + string_of_int(i+33), "J" + string_of_int(i+33), j_seq[i], "", NULL};
293 294 295 296 297 298 299 300 301 302 303 304 305
    V.add(v);
    J.add(j);
  }
  Germline germline("Test", 'T', V, Fasta(), J, 0, 30);
  germline.new_index("####");


  TAP_TEST(germline.index->getIndexLoad() == .75, TEST_GET_INDEX_LOAD, "");


  Sequence seq = {"to_segment", "to_segment", "TATCG", "", NULL};
  KmerSegmenter kseg(seq, &germline);

306 307 308 309 310
  KmerAffectAnalyser *kaa = kseg.getKmerAffectAnalyser();

  TAP_TEST(kaa->getProbabilityAtLeastOrAbove(3) == 0, TEST_PROBABILITY_SEGMENTATION, "");
  TAP_TEST(kaa->getProbabilityAtLeastOrAbove(2) == .75 * .75, TEST_PROBABILITY_SEGMENTATION, "");
  TAP_TEST(kaa->getProbabilityAtLeastOrAbove(1) == .75 * 2 * .25 + kaa->getProbabilityAtLeastOrAbove(2),
311
            TEST_PROBABILITY_SEGMENTATION, "");
312
  TAP_TEST(kaa->getProbabilityAtLeastOrAbove(0) == 1, TEST_PROBABILITY_SEGMENTATION, "");
313 314
}

Mikael Salson's avatar
Mikael Salson committed
315 316 317 318
void testSegment() {
  testSegmentOverlap();
  testSegmentationCause();
  testExtractor();
319
  testProbability();
320
  testFineSegment();
Mikael Salson's avatar
Mikael Salson committed
321
}