Medical Imaging Interaction Toolkit  2016.11.0
Medical Imaging Interaction Toolkit
mitkTumorInvasionClassification.cpp
Go to the documentation of this file.
1 /*===================================================================
2 
3 The Medical Imaging Interaction Toolkit (MITK)
4 
5 Copyright (c) German Cancer Research Center,
6 Division of Medical and Biological Informatics.
7 All rights reserved.
8 
9 This software is distributed WITHOUT ANY WARRANTY; without
10 even the implied warranty of MERCHANTABILITY or FITNESS FOR
11 A PARTICULAR PURPOSE.
12 
13 See LICENSE.txt or http://www.mitk.org for details.
14 
15 ===================================================================*/
16 
17 #include "mitkTumorInvasionClassification.h"
18 
19 #include "mitkCLUtil.h"
20 #include "mitkCollectionDilatation.h"
21 #include "mitkCollectionGrayOpening.h"
22 #include "mitkDataCollectionImageIterator.h"
23 #include "mitkDataCollectionUtilities.h"
24 #include "mitkIOUtil.h"
25 #include <fstream>
26 // To initialize random number generator
27 #include <time.h>
28 
29 static void EnsureDataImageInCollection(mitk::DataCollection::Pointer collection,
30  std::string origin,
31  std::string target)
32 {
33  typedef itk::Image<double, 3> FeatureImage;
34 
35  if (collection->HasElement(origin))
36  {
37  mitk::Image::Pointer originImage = dynamic_cast<mitk::Image *>(collection->GetMitkImage(origin).GetPointer());
38  FeatureImage::Pointer itkOriginImage = FeatureImage::New();
39  mitk::CastToItkImage(originImage, itkOriginImage);
40 
41  if (!collection->HasElement(target) && itkOriginImage.IsNotNull())
42  {
43  FeatureImage::Pointer image = FeatureImage::New();
44  image->SetRegions(itkOriginImage->GetLargestPossibleRegion());
45  image->SetSpacing(itkOriginImage->GetSpacing());
46  image->SetOrigin(itkOriginImage->GetOrigin());
47  image->SetDirection(itkOriginImage->GetDirection());
48  image->Allocate();
49 
50  collection->AddData(dynamic_cast<itk::DataObject *>(image.GetPointer()), target, "");
51  }
52  }
53  for (std::size_t i = 0; i < collection->Size(); ++i)
54  {
55  mitk::DataCollection *newCol = dynamic_cast<mitk::DataCollection *>(collection->GetData(i).GetPointer());
56  if (newCol != 0)
57  {
58  EnsureDataImageInCollection(newCol, origin, target);
59  }
60  }
61 }
62 
63 mitk::TumorInvasionClassification::TumorInvasionClassification()
64  : m_TargetID("TARGET"),
65  m_TumorID("GTV"),
66  m_MaskID("BRAINMASK"),
67  m_ResultID("RESULT"),
68  m_Randomize(false),
69  m_WeightSamples(false)
70 {
71 }
72 
73 void mitk::TumorInvasionClassification::SelectTrainingSamples(mitk::DataCollection *collection, unsigned int mode)
74 {
75  srand(time(NULL));
76  MITK_INFO << "LearnProgressionFeatures: Selecting training voxels.";
77  switch (mode)
78  {
79  case 0:
80  {
81  MITK_INFO << " Selection Mode " << mode << " use all tumor voxels, healthy: 50% vicinity / 50% far away";
82 
83  CollectionDilation::DilateBinaryByName(collection, m_TargetID, 2, 0, "EXCLUDE");
84  CollectionDilation::ErodeBinaryByName(collection, m_TargetID, 1, 0, "ERODE");
85  CollectionDilation::DilateBinaryByName(collection, m_TargetID, m_TargetDil2D, m_TargetDil3D, "TRAIN");
86  DataCollectionImageIterator<unsigned char, 3> gtvIter(collection, m_TumorID);
87  DataCollectionImageIterator<unsigned char, 3> brainMaskIter(collection, m_MaskID);
88  DataCollectionImageIterator<unsigned char, 3> targetIter(collection, m_TargetID);
89  DataCollectionImageIterator<unsigned char, 3> targetDil(collection, m_TargetID + "TRAIN");
90 
91  // Count Healthy/ Tumor voxels
92  // i.o. to decide how many of each are taken
93  unsigned int totalTumor = 0;
94  unsigned int totalHealthyClose = 0;
95  unsigned int totalHealthy = 0;
96 
97  while (!brainMaskIter.IsAtEnd())
98  {
99  if (brainMaskIter.GetVoxel() != 0)
100  {
101  if (targetIter.GetVoxel() == 1 && gtvIter.GetVoxel() == 0)
102  ++totalTumor;
103 
104  if (targetIter.GetVoxel() == 0 && targetDil.GetVoxel() == 1 && gtvIter.GetVoxel() == 0)
105  ++totalHealthyClose;
106 
107  if (targetIter.GetVoxel() == 0 && gtvIter.GetVoxel() == 0 && targetDil.GetVoxel() == 0)
108  ++totalHealthy; // healthy but not close
109  }
110  ++brainMaskIter;
111  ++targetIter;
112  ++targetDil;
113  ++gtvIter;
114  }
115  brainMaskIter.ToBegin();
116  targetIter.ToBegin();
117  targetDil.ToBegin();
118  gtvIter.ToBegin();
119  // Total of healthy samples that is to be collected
120  unsigned int targetHealthy = totalTumor * m_ClassRatio;
121  // Determines which portion of the healthy samples is drawn from the immediate vicinity of the newly grown tumor
122  ScalarType ratioClose = .5;
123 
124  // Compute probabilities thresholds for choosing a close healthy voxel / any healthy voxel
125  ScalarType thHealthyClose = std::min(1.0, (targetHealthy * ratioClose) / totalHealthyClose);
126  ScalarType thHealthyAny = std::min(1.0, (targetHealthy * (1.0 - ratioClose)) / totalHealthy);
127  // Some stats
128  {
129  MITK_INFO << "Total # Tumor Voxels " << totalTumor;
130  MITK_INFO << "Total # Healthy Voxels" << totalHealthy;
131  MITK_INFO << "Target Ratio " << m_ClassRatio;
132  MITK_INFO << "Ratio for healthy close: " << thHealthyClose;
133  MITK_INFO << "Ratio for healthy any other: " << thHealthyAny;
134  }
135 
136  // DEBUG count occurances to compare
137  unsigned int tumor = 0;
138  unsigned int healthy = 0;
139  while (!brainMaskIter.IsAtEnd())
140  {
141  if (brainMaskIter.GetVoxel() != 0)
142  {
143  if (targetIter.GetVoxel() == 1)
144  { // choose tumor voxels for training
145  if (gtvIter.GetVoxel() == 0) // tumor always
146  {
147  targetIter.SetVoxel(brainMaskIter.GetVoxel() + targetIter.GetVoxel());
148  ++tumor;
149  }
150  else
151  targetIter.SetVoxel(0);
152  }
153  else
154  { // choose healty tissue voxels for training
155  ScalarType rndVal = (float)rand() / (float)(RAND_MAX); //(0..1)
156  if (gtvIter.GetVoxel() == 0 && ((targetDil.GetVoxel() == 1 && rndVal <= thHealthyClose)))
157  {
158  targetIter.SetVoxel(brainMaskIter.GetVoxel() + targetIter.GetVoxel());
159  ++healthy;
160  }
161  else if (((targetDil.GetVoxel() == 0 && rndVal <= thHealthyAny)))
162  {
163  targetIter.SetVoxel(brainMaskIter.GetVoxel() + targetIter.GetVoxel());
164  ++healthy;
165  }
166  else
167  targetIter.SetVoxel(0);
168  }
169  }
170  else
171  targetIter.SetVoxel(0);
172 
173  ++brainMaskIter;
174  ++targetIter;
175  ++gtvIter;
176  ++targetDil;
177  }
178  MITK_INFO << "Training with Samples #Tumor " << tumor << " / healthy # " << healthy;
179  }
180  break;
181  case 2:
182  {
183  MITK_INFO << " Selection Mode " << mode << " Weighted with ratio";
184 
185  EnsureDataImageInCollection(collection, m_TumorID, "WEIGHTS");
186 
187  CollectionDilation::DilateBinaryByName(collection, m_TargetID, 1, 0, "EXCLUDE");
188  CollectionDilation::ErodeBinaryByName(collection, m_TargetID, 1, 0, "ERODE");
189  CollectionDilation::DilateBinaryByName(collection, m_TargetID + "EXCLUDE", m_TargetDil2D, m_TargetDil3D, "TRAIN");
190  DataCollectionImageIterator<unsigned char, 3> gtvIter(collection, m_TumorID);
191  DataCollectionImageIterator<unsigned char, 3> brainMaskIter(collection, m_MaskID);
192  DataCollectionImageIterator<unsigned char, 3> targetIter(collection, m_TargetID);
193  DataCollectionImageIterator<unsigned char, 3> targetDil(collection, m_TargetID + "EXCLUDETRAIN");
194  DataCollectionImageIterator<unsigned char, 3> excludeTumorIter(collection, m_TargetID + "ERODE");
195  DataCollectionImageIterator<unsigned char, 3> excludeHealthyIter(collection, m_TargetID + "EXCLUDE");
196 
197  DataCollectionImageIterator<double, 3> weightsIter(collection, "WEIGHTS");
198 
199  // Count Healthy/ Tumor voxels
200  // i.o. to decide how many of each are taken
201  double totalTumor = 0;
202  double totalHealthyClose = 0;
203  double totalHealthyNonClose = 0;
204  double totalHealthy = 0;
205 
206  while (!brainMaskIter.IsAtEnd())
207  {
208  if (brainMaskIter.GetVoxel() != 0)
209  {
210  if (targetIter.GetVoxel() == 1 && gtvIter.GetVoxel() == 0 && excludeTumorIter.GetVoxel() == 1)
211  ++totalTumor;
212 
213  if (excludeHealthyIter.GetVoxel() == 0 && targetDil.GetVoxel() == 1 && gtvIter.GetVoxel() == 0)
214  ++totalHealthyClose;
215 
216  if (excludeHealthyIter.GetVoxel() == 0 && gtvIter.GetVoxel() == 0 && targetDil.GetVoxel() == 0)
217  ++totalHealthyNonClose; // healthy but not close
218 
219  if (excludeHealthyIter.GetVoxel() == 0 && gtvIter.GetVoxel() == 0 && targetIter.GetVoxel() == 0)
220  ++totalHealthy; // healthy
221  }
222  ++brainMaskIter;
223  ++targetIter;
224  ++targetDil;
225  ++gtvIter;
226  ++excludeHealthyIter;
227  ++excludeTumorIter;
228  }
229  brainMaskIter.ToBegin();
230  targetIter.ToBegin();
231  targetDil.ToBegin();
232  gtvIter.ToBegin();
233  excludeHealthyIter.ToBegin();
234  excludeTumorIter.ToBegin();
235 
236  // Compute probabilities thresholds for choosing a close healthy voxel / any healthy voxel
237  ScalarType wHealthyClose = 10000.0 / totalHealthyClose;
238  ScalarType wHealthyAny = 5000.0 / totalHealthyNonClose;
239  ScalarType wTumor = (m_ClassRatio * (10000.0 + 5000.0)) / totalTumor;
240 
241  // DEBUG count occurances to compare
242  double potentialClose = 0;
243  double selectedClose = 0;
244  double tumor = 0;
245  double healthy = 0;
246  while (!brainMaskIter.IsAtEnd())
247  {
248  weightsIter.SetVoxel(0);
249  if (brainMaskIter.GetVoxel() != 0)
250  {
251  if (targetIter.GetVoxel() == 1)
252  { // choose tumor voxels for training
253  if (gtvIter.GetVoxel() == 0 && excludeTumorIter.GetVoxel() == 1) // tumor always
254  {
255  targetIter.SetVoxel(brainMaskIter.GetVoxel() + targetIter.GetVoxel());
256  weightsIter.SetVoxel(wTumor);
257  ++tumor;
258  }
259  else if (gtvIter.GetVoxel() == 0)
260  {
261  weightsIter.SetVoxel(0); //.1);
262  targetIter.SetVoxel(0); // 2);
263  }
264  }
265  else
266  { // choose healty tissue voxels for training
267  if (gtvIter.GetVoxel() == 0 && ((excludeHealthyIter.GetVoxel() == 0 && targetDil.GetVoxel() == 1)))
268  {
269  targetIter.SetVoxel(brainMaskIter.GetVoxel() + targetIter.GetVoxel());
270  weightsIter.SetVoxel(wHealthyClose);
271  healthy += wHealthyClose;
272  }
273  else if (((targetDil.GetVoxel() == 0 && excludeHealthyIter.GetVoxel() == 0)))
274  {
275  targetIter.SetVoxel(brainMaskIter.GetVoxel() + targetIter.GetVoxel());
276  weightsIter.SetVoxel(wHealthyAny);
277  healthy += wHealthyAny;
278  }
279  else if ((gtvIter.GetVoxel() == 0 && excludeHealthyIter.GetVoxel() == 1))
280  {
281  targetIter.SetVoxel(0); // brainMaskIter.GetVoxel()+targetIter.GetVoxel());
282  weightsIter.SetVoxel(0); //.1);
283  }
284  else
285  {
286  targetIter.SetVoxel(0);
287  weightsIter.SetVoxel(0);
288  }
289  }
290  }
291  else
292  {
293  targetIter.SetVoxel(0);
294  weightsIter.SetVoxel(0);
295  }
296 
297  if (gtvIter.GetVoxel() != 0)
298  targetIter.SetVoxel(0);
299 
300  ++brainMaskIter;
301  ++targetIter;
302  ++gtvIter;
303  ++targetDil;
304  ++excludeHealthyIter;
305  ++excludeTumorIter;
306  ++weightsIter;
307  }
308  MITK_INFO << "Training with Samples #Tumor " << tumor << " / healthy # " << healthy;
309  MITK_INFO << "Potential Close" << potentialClose;
310  MITK_INFO << "Selected Close" << selectedClose;
311  }
312  break;
313  default:
314  {
315  MITK_INFO << " Selection Mode " << mode
316  << " Exclude voxels in border regions, healthy: 50% vicinity / 50% far away";
317 
318  // weights
319  ScalarType tumorWeight = 1;
320  // ScalarType unsureRegion = .25;
321  ScalarType healthyTissue = 1;
322 
323  EnsureDataImageInCollection(collection, m_TumorID, "WEIGHTS");
324 
325  CollectionDilation::DilateBinaryByName(collection, m_TargetID, 1, 0, "EXCLUDE");
326  CollectionDilation::ErodeBinaryByName(collection, m_TargetID, 1, 0, "ERODE");
327  CollectionDilation::DilateBinaryByName(collection, m_TargetID + "EXCLUDE", m_TargetDil2D, m_TargetDil3D, "TRAIN");
328  DataCollectionImageIterator<unsigned char, 3> gtvIter(collection, m_TumorID);
329  DataCollectionImageIterator<unsigned char, 3> brainMaskIter(collection, m_MaskID);
330  DataCollectionImageIterator<unsigned char, 3> targetIter(collection, m_TargetID);
331  DataCollectionImageIterator<unsigned char, 3> targetDil(collection, m_TargetID + "EXCLUDETRAIN");
332  DataCollectionImageIterator<unsigned char, 3> excludeTumorIter(collection, m_TargetID + "ERODE");
333  DataCollectionImageIterator<unsigned char, 3> excludeHealthyIter(collection, m_TargetID + "EXCLUDE");
334 
335  DataCollectionImageIterator<double, 3> weightsIter(collection, "WEIGHTS");
336 
337  // Count Healthy/ Tumor voxels
338  // i.o. to decide how many of each are taken
339  double totalTumor = 0;
340  double totalHealthyClose = 0;
341  double totalHealthyNonClose = 0;
342  double totalHealthy = 0;
343 
344  while (!brainMaskIter.IsAtEnd())
345  {
346  if (brainMaskIter.GetVoxel() != 0)
347  {
348  if (targetIter.GetVoxel() == 1 && gtvIter.GetVoxel() == 0 && excludeTumorIter.GetVoxel() == 1)
349  ++totalTumor;
350 
351  if (excludeHealthyIter.GetVoxel() == 0 && targetDil.GetVoxel() == 1 && gtvIter.GetVoxel() == 0)
352  ++totalHealthyClose;
353 
354  if (excludeHealthyIter.GetVoxel() == 0 && gtvIter.GetVoxel() == 0 && targetDil.GetVoxel() == 0)
355  ++totalHealthyNonClose; // healthy but not close
356 
357  if (excludeHealthyIter.GetVoxel() == 0 && gtvIter.GetVoxel() == 0 && targetIter.GetVoxel() == 0)
358  ++totalHealthy; // healthy
359  }
360  ++brainMaskIter;
361  ++targetIter;
362  ++targetDil;
363  ++gtvIter;
364  ++excludeHealthyIter;
365  ++excludeTumorIter;
366  }
367  brainMaskIter.ToBegin();
368  targetIter.ToBegin();
369  targetDil.ToBegin();
370  gtvIter.ToBegin();
371  excludeHealthyIter.ToBegin();
372  excludeTumorIter.ToBegin();
373 
374  // Total of healthy samples that is to be collected
375  unsigned int targetHealthy = (tumorWeight / healthyTissue) * totalTumor * m_ClassRatio;
376  // Determines which portion of the healthy samples is drawn from the immediate vicinity of the newly grown tumor
377  ScalarType ratioClose = .5;
378 
379  // Compute probabilities thresholds for choosing a close healthy voxel / any healthy voxel
380  ScalarType thHealthyClose = std::min(1.0, ((double)targetHealthy * ratioClose) / totalHealthyClose);
381  ScalarType thHealthyAny = std::min(1.0, ((double)targetHealthy * (1.0 - ratioClose)) / totalHealthyNonClose);
382  // Some stats
383  {
384  MITK_INFO << "Total Tumor " << totalTumor;
385  MITK_INFO << "Total healthy " << totalHealthyNonClose;
386  MITK_INFO << "Total healthy close " << totalHealthyClose;
387  MITK_INFO << "Total healthy non-close " << totalHealthyNonClose;
388  MITK_INFO << "Target Ratio " << m_ClassRatio;
389 
390  MITK_INFO << "Target Healthy " << targetHealthy;
391  MITK_INFO << "Probabilty close " << thHealthyClose;
392  MITK_INFO << "Probabilty any other " << thHealthyAny;
393  }
394 
395  // DEBUG count occurances to compare
396  double potentialClose = 0;
397  double selectedClose = 0;
398  double tumor = 0;
399  double healthy = 0;
400  while (!brainMaskIter.IsAtEnd())
401  {
402  weightsIter.SetVoxel(0);
403  if (brainMaskIter.GetVoxel() != 0)
404  {
405  if (targetIter.GetVoxel() == 1)
406  { // choose tumor voxels for training
407  if (gtvIter.GetVoxel() == 0 && excludeTumorIter.GetVoxel() == 1) // tumor always
408  {
409  targetIter.SetVoxel(brainMaskIter.GetVoxel() + targetIter.GetVoxel());
410  weightsIter.SetVoxel(tumorWeight);
411  ++tumor;
412  }
413  else if (gtvIter.GetVoxel() == 0)
414  {
415  weightsIter.SetVoxel(0); //.1);
416  targetIter.SetVoxel(0); // 2);
417  }
418  }
419  else
420  { // choose healty tissue voxels for training
421  if (gtvIter.GetVoxel() == 0 && ((excludeHealthyIter.GetVoxel() == 0 && targetDil.GetVoxel() == 1)))
422  {
423  if (!m_WeightSamples && ((float)rand() / (float)(RAND_MAX) < thHealthyClose))
424  {
425  targetIter.SetVoxel(brainMaskIter.GetVoxel() + targetIter.GetVoxel());
426  weightsIter.SetVoxel(1);
427  ++healthy;
428  }
429  else
430  {
431  targetIter.SetVoxel(brainMaskIter.GetVoxel() + targetIter.GetVoxel());
432  weightsIter.SetVoxel(thHealthyClose);
433  healthy += thHealthyClose;
434  }
435  }
436  else if (((targetDil.GetVoxel() == 0 && excludeHealthyIter.GetVoxel() == 0)))
437  {
438  if (!m_WeightSamples && ((float)rand() / (float)(RAND_MAX) < thHealthyAny))
439  {
440  targetIter.SetVoxel(brainMaskIter.GetVoxel() + targetIter.GetVoxel());
441  weightsIter.SetVoxel(1);
442  ++healthy;
443  }
444  else
445  {
446  targetIter.SetVoxel(brainMaskIter.GetVoxel() + targetIter.GetVoxel());
447  weightsIter.SetVoxel(thHealthyAny);
448  healthy += thHealthyAny;
449  }
450  }
451  else if ((gtvIter.GetVoxel() == 0 && excludeHealthyIter.GetVoxel() == 1))
452  {
453  targetIter.SetVoxel(0); // brainMaskIter.GetVoxel()+targetIter.GetVoxel());
454  weightsIter.SetVoxel(0); //.1);
455  }
456  else
457  {
458  targetIter.SetVoxel(0);
459  weightsIter.SetVoxel(0);
460  }
461 
462  if (gtvIter.GetVoxel() == 0 && ((excludeHealthyIter.GetVoxel() == 0 && targetDil.GetVoxel() == 1)))
463  {
464  potentialClose++;
465 
466  if ((float)rand() / (float)(RAND_MAX) < thHealthyClose)
467  selectedClose++;
468  }
469  }
470  }
471  else
472  {
473  targetIter.SetVoxel(0);
474  weightsIter.SetVoxel(0);
475  }
476 
477  if (gtvIter.GetVoxel() != 0)
478  targetIter.SetVoxel(0);
479 
480  ++brainMaskIter;
481  ++targetIter;
482  ++gtvIter;
483  ++targetDil;
484  ++excludeHealthyIter;
485  ++excludeTumorIter;
486  ++weightsIter;
487  }
488  MITK_INFO << "Training with Samples #Tumor " << tumor << " / healthy # " << healthy;
489  MITK_INFO << "Potential Close" << potentialClose;
490  MITK_INFO << "Selected Close" << selectedClose;
491  }
492  break;
493  }
494 }
495 
496 void mitk::TumorInvasionClassification::PrepareResponseSamples(mitk::DataCollection *collection)
497 {
498  srand(time(NULL));
499  MITK_INFO << "PrepareResponseSamples: Selecting training voxels.";
500 
501  EnsureDataImageInCollection(collection, m_TumorID, "WEIGHTS");
502 
503  CollectionDilation::DilateBinaryByName(collection, m_TargetID, 10, 4, "TRAIN");
504  DataCollectionImageIterator<unsigned char, 3> gtvIter(collection, m_TumorID);
505  DataCollectionImageIterator<unsigned char, 3> brainMaskIter(collection, m_MaskID);
506  DataCollectionImageIterator<unsigned char, 3> targetIter(collection, m_TargetID);
507  DataCollectionImageIterator<unsigned char, 3> dilIter(collection, m_TargetID + "TRAIN");
508 
509  // Count Non-Involved/Responsive/Rezidiv Voxels
510 
511  double nonInvolved = 0;
512  double responsive = 0;
513  double rezidiv = 0;
514 
515  while (!brainMaskIter.IsAtEnd())
516  {
517  if (brainMaskIter.GetVoxel() != 0 && dilIter.GetVoxel() != 0)
518  {
519  if (targetIter.GetVoxel() == 0 && gtvIter.GetVoxel() == 0)
520  {
521  targetIter.SetVoxel(1);
522  ++nonInvolved;
523  }
524 
525  if (targetIter.GetVoxel() == 0 && gtvIter.GetVoxel() == 1)
526  {
527  targetIter.SetVoxel(3);
528  ++responsive;
529  }
530 
531  if (targetIter.GetVoxel() == 1 && gtvIter.GetVoxel() == 0)
532  {
533  targetIter.SetVoxel(2);
534  ++rezidiv;
535  }
536  }
537  else
538  targetIter.SetVoxel(0);
539 
540  ++brainMaskIter;
541  ++targetIter;
542  ++gtvIter;
543  ++dilIter;
544  }
545  brainMaskIter.ToBegin();
546  targetIter.ToBegin();
547  gtvIter.ToBegin();
548  dilIter.ToBegin();
549 
550  // weight for groups
551 
552  double wNonInvolved = 10000.0 / nonInvolved;
553  double wResponsive = 10000.0 / responsive;
554  double wRezidiv = 10000.0 / rezidiv;
555 
556  std::cout << "Weights are " << wNonInvolved << "/ " << wResponsive << " / " << wRezidiv << std::endl;
557 
558  // Assign weight for each voxel
559  DataCollectionImageIterator<double, 3> weightsIter(collection, "WEIGHTS");
560 
561  while (!brainMaskIter.IsAtEnd())
562  {
563  if (brainMaskIter.GetVoxel() != 0 && dilIter.GetVoxel() != 0)
564  {
565  if (targetIter.GetVoxel() == 1)
566  weightsIter.SetVoxel(wNonInvolved);
567 
568  if (targetIter.GetVoxel() == 2)
569  weightsIter.SetVoxel(wRezidiv);
570 
571  if (targetIter.GetVoxel() == 3)
572  weightsIter.SetVoxel(wResponsive);
573  }
574 
575  ++dilIter;
576  ++brainMaskIter;
577  ++targetIter;
578  ++weightsIter;
579  }
580 }
581 
582 void mitk::TumorInvasionClassification::LearnProgressionFeatures(mitk::DataCollection *collection,
583  std::vector<std::string> modalitiesList,
584  size_t forestSize,
585  size_t treeDepth)
586 {
587  auto trainDataX = mitk::DCUtilities::DC3dDToMatrixXd(collection, modalitiesList, m_TargetID);
588  auto trainDataY = mitk::DCUtilities::DC3dDToMatrixXi(collection, m_TargetID, m_TargetID);
589 
590  if (treeDepth != 0)
591  m_Forest.SetMaximumTreeDepth(treeDepth);
592 
593  m_Forest.SetTreeCount(forestSize);
594 
595  if (m_WeightSamples)
596  {
597  auto trainDataW = mitk::DCUtilities::DC3dDToMatrixXd(collection, "WEIGHTS", m_TargetID);
598  m_Forest.SetPointWiseWeight(trainDataW);
599  m_Forest.UsePointWiseWeight(true);
600  }
601 
602  MITK_INFO << "Start Training:";
603  try
604  {
605  m_Forest.Train(trainDataX, trainDataY);
606  }
607  catch (std::exception &e)
608  {
609  MITK_INFO << "Exception while training forest: " << e.what();
610  }
611  std::cout << "Training finished" << std::endl;
612 }
613 
614 void mitk::TumorInvasionClassification::PredictInvasion(mitk::DataCollection *collection,
615  std::vector<std::string> modalitiesList)
616 {
617  if (collection != NULL)
618  {
619  MITK_INFO << "Evaluating Forest";
620 
621  auto testDataX = mitk::DCUtilities::DC3dDToMatrixXd(collection, modalitiesList, m_MaskID);
622  auto testDataNewY = m_Forest.Predict(testDataX);
623  mitk::DCUtilities::MatrixToDC3d(testDataNewY, collection, m_ResultID, m_MaskID);
624 
625  Eigen::MatrixXd Probs = m_Forest.GetPointWiseProbabilities();
626 
627  Eigen::MatrixXd prob0 = Probs.col(0);
628  Eigen::MatrixXd prob1 = Probs.col(1);
629 
630  mitk::DCUtilities::MatrixToDC3d(prob0, collection, "prob0", m_MaskID);
631  mitk::DCUtilities::MatrixToDC3d(prob1, collection, "prob1", m_MaskID);
632  if (Probs.cols() >= 3)
633  {
634  Eigen::MatrixXd prob2 = Probs.col(2);
635  mitk::DCUtilities::MatrixToDC3d(prob2, collection, "prob2", m_MaskID);
636  }
637  }
638  else
639  MITK_ERROR << "TumorInvasionClassification::PredictInvasion - provided collection is NULL.";
640 }
641 
642 void mitk::TumorInvasionClassification::SanitizeResults(mitk::DataCollection *collection, std::string resultID)
643 {
644  CollectionGrayOpening::PerformGrayOpening(collection, resultID, "OPEN");
645 }
646 
647 void mitk::TumorInvasionClassification::DescriptionToLogFile(mitk::DataCollection *collection, std::string outputFile)
648 {
649  std::ofstream fileStream;
650  fileStream.open(outputFile.c_str(), std::ios_base::app);
651 
652  fileStream << "Configuration:" << std::endl;
653  for (size_t patient = 0; patient < collection->Size(); ++patient)
654  {
655  DataCollection *dataPatient = dynamic_cast<DataCollection *>(collection->GetData(patient).GetPointer());
656  fileStream << "Name : " << dataPatient->GetName() << " with time steps:" << std::endl;
657  for (size_t timeStep = 0; timeStep < dataPatient->Size(); ++timeStep)
658  {
659  fileStream << dataPatient->IndexToName(timeStep) << "-" << std::endl;
660  DataCollection *dataTS = dynamic_cast<DataCollection *>(dataPatient->GetData(timeStep).GetPointer());
661  for (size_t ts = 0; ts < dataTS->Size(); ++ts)
662  fileStream << dataTS->IndexToName(ts) << "-" << std::endl;
663  }
664  }
665 }
666 
667 void mitk::TumorInvasionClassification::SaveRandomForest(std::string filename)
668 {
669  mitk::IOUtil::Save(&m_Forest, filename);
670 }
671 
672 void mitk::TumorInvasionClassification::LoadRandomForest(std::string /*filename*/)
673 {
674  // TBD
675  MITK_ERROR << "not yet implemented ...";
676 }
mitk::CollectionDilation::ErodeBinaryByName
static void ErodeBinaryByName(DataCollection *dataCollection, std::string name, unsigned int xy=5, unsigned int z=0, std::string suffix="ERODE")
ErodeBinaryByName - Erode all occurances of a modality within a mitk::DataCollection.
Definition: mitkCollectionDilatation.cpp:76
mitk::IOUtil::Save
static void Save(const mitk::BaseData *data, const std::string &path)
Save a mitk::BaseData instance.
Definition: mitkIOUtil.cpp:824
FeatureImage
itk::Image< mitk::ScalarType, 3 > FeatureImage
Definition: HistogramAdaption.cpp:26
mitk::Pointer
itk::SmartPointer< Self > Pointer
Definition: mitkRenderingManager.h:389
MITK_INFO
#define MITK_INFO
Definition: mitkLogMacros.h:22
mitk::TumorInvasionClassification::SaveRandomForest
void SaveRandomForest(std::string filename)
SaveRandomForest - Saves a trained random forest.
Definition: mitkTumorInvasionClassification.cpp:667
MITK_ERROR
#define MITK_ERROR
Definition: mitkLogMacros.h:24
mitk::ScalarType
double ScalarType
Definition: mitkNumericConstants.h:24
mitk::DataCollection::Size
vcl_size_t Size() const
Size - number of data items in collection.
Definition: mitkDataCollection.cpp:116
mitk::DataCollection::GetName
std::string GetName() const
Definition: mitkDataCollection.cpp:61
mitk::TumorInvasionClassification::DescriptionToLogFile
void DescriptionToLogFile(DataCollection *collection, std::string outputFile)
DescriptionToLogFile - Append to log file: patient ids and time steps contained in this data collecti...
Definition: mitkTumorInvasionClassification.cpp:647
mitk::TumorInvasionClassification::SelectTrainingSamples
void SelectTrainingSamples(DataCollection *collection, unsigned int mode=0)
SelectTrainingSamples.
Definition: mitkTumorInvasionClassification.cpp:73
mitk::TumorInvasionClassification::PrepareResponseSamples
void PrepareResponseSamples(DataCollection *collection)
PrepareResponseSamples.
Definition: mitkTumorInvasionClassification.cpp:496
mitk::TumorInvasionClassification::PredictInvasion
void PredictInvasion(DataCollection *collection, std::vector< std::string > modalitiesList)
PredictGrowth - Classify voxels into remaining healthy / turning into tumor.
Definition: mitkTumorInvasionClassification.cpp:614
mitk::DataCollection::IndexToName
std::string IndexToName(vcl_size_t index) const
IndexToName - Get name from index.
Definition: mitkDataCollection.cpp:90
mitk::CollectionGrayOpening::PerformGrayOpening
static void PerformGrayOpening(mitk::DataCollection *dataCollection, std::string name, std::string suffix)
PerformGrayOpening - Opening operation on a specific modality type wihtin the DataCollection. Creates a new item.
Definition: mitkCollectionGrayOpening.cpp:26
mitk::DCUtilities::MatrixToDC3d
static void MatrixToDC3d(const Eigen::MatrixXd &matrix, mitk::DataCollection::Pointer dc, const std::vector< std::string > &names, std::string mask)
Definition: mitkDataCollectionUtilities.cpp:125
filename
static const std::string filename
Definition: mitkTinyXMLTest.cpp:28
mitk::DCUtilities::DC3dDToMatrixXi
static Eigen::MatrixXi DC3dDToMatrixXi(mitk::DataCollection::Pointer dc, std::string name, std::string mask)
Definition: mitkDataCollectionUtilities.cpp:80
mitk::Image
Image class for storing images.
Definition: mitkImage.h:76
mitk::DCUtilities::DC3dDToMatrixXd
static Eigen::MatrixXd DC3dDToMatrixXd(mitk::DataCollection::Pointer dc, std::string names, std::string mask)
Definition: mitkDataCollectionUtilities.cpp:36
min
static T min(T x, T y)
Definition: svm.cpp:67
mitk::TumorInvasionClassification::LearnProgressionFeatures
void LearnProgressionFeatures(DataCollection *collection, std::vector< std::string > modalitiesList, vcl_size_t forestSize=300, vcl_size_t treeDepth=10)
LearnProgressionFeatures.
Definition: mitkTumorInvasionClassification.cpp:582
mitk::CastToItkImage
void MITKCORE_EXPORT CastToItkImage(const mitk::Image *mitkImage, itk::SmartPointer< ItkOutputImageType > &itkOutputImage)
Cast an mitk::Image to an itk::Image with a specific type.
EnsureDataImageInCollection
static void EnsureDataImageInCollection(mitk::DataCollection::Pointer collection, std::string origin, std::string target)
Definition: mitkTumorInvasionClassification.cpp:29
mitk::DataCollection::GetData
itk::DataObject::Pointer GetData(vcl_size_t index)
GetData Get original data by index.
Definition: mitkDataCollection.cpp:141
mitk::TumorInvasionClassification::SanitizeResults
void SanitizeResults(DataCollection *collection, std::string resultID="RESULT")
SanitizeResults - Performs an Opening Operation on tha data to remove isolated misclassifications.
Definition: mitkTumorInvasionClassification.cpp:642
mitk::TumorInvasionClassification::LoadRandomForest
void LoadRandomForest(std::string)
LoadRandomForest - loads a random forest to be used for classification.
Definition: mitkTumorInvasionClassification.cpp:672
mitk::CollectionDilation::DilateBinaryByName
static void DilateBinaryByName(DataCollection *dataCollection, std::string name, unsigned int xy=5, unsigned int z=0, std::string suffix="DILATE")
DilateBinaryByName - Dilate all occurances of a modality within a mitk::DataCollection.
Definition: mitkCollectionDilatation.cpp:27
mitk::New
static itkEventMacro(BoundingShapeInteractionEvent, itk::AnyEvent) class MITKBOUNDINGSHAPE_EXPORT BoundingShapeInteractor Pointer New()
Basic interaction methods for mitk::GeometryData.