There are several reasons why we rarely see 3D agent based models of tumor growth, with I think the two most important ones being: 1) 3D simulations are way more computationally expensive than 2D simulations; 2) it’s harder to visualize results of 3D simulations to present all necessary information. The other frequently used argument is that there is no need for 3D, because the crucial phenomena can be observed as well in 2D. However, when one considers e.g. stochastic models there might be a huge difference in the system behavior between 2D and 3D case (see very nice post Why is that a 2D random walk is recurrent, while a 3D random walk is transient?).
Possible improvements in the computational speed are highly dependent on the given model setting and typically require sophisticated techniques. Problem with nice visualizations, however, just requires proper tools. Today, I will show my attempts to make a good-looking 3D tumor visualizations, which I have taken over the last several years.
Several years ago when working with @heikman on the project about the impact of cellular senescence I was using The Persistence of Vision Raytracer (POV-Ray) in order to obtain high quality 3D tumor visualization. After simulating the tumor up to a given cell number I was parsing the information about cells locations and their basic properties into a POV-Ray specific format. In my approach, each cell was separate actor in the scene with possible different surface/texture/color properties. As it can be expected renderings with raytracing software are very good (see Fig.1 below), but the whole procedure to prepare input files to POV-Ray, to set-up properly camera and light was really time-consuming. In addition, rendering takes a long time, so there is no chance for “live” rendering during the simulation, which would allow to observe simulation behavior.
Figure 1. Exemplary 3D tumor renderings obtained using POV-Ray software. Bottom-left picture shows tumor consisted of 9.4 million of cells.
I really wanted to have a tool that would allow me to seamlessly simulate and visualize 3D tumor at the same time. I quickly realized that rendering each cell as a separate object is not a good approach, when you simulate more than 1 million cells. Hence, I started looking around for an algorithm that would allow to extract tumor surface from the cloud of cells. I quickly found well known Marching Cubes algorithm that does the trick very well. I’ve implemented the algorithm and visualized resulting surface using Visualization Library in C++. I was quite happy with the results (see Fig. 2 below), however, at some point code became so complicated and so model specific that I assumed that it doesn’t have any future.
Figure 2. “Live” 3D tumor visualization during the simulation. Code developed in C++ using Visualization Library.
Quite recently I started using VTK library from within C++ (you can also use it form within Python) and I quickly realized that it might be the ultimate tool for 3D visualization. Now I will show how easy it is to produce great visualizations that can be used in “live” renderings during 3D simulations using VTK library.
In the following example I will use glyph mapper that basically takes my cloud of points (cells), creates cubes around them, and extracts the surface. First, let us include all necessary headers and define necessary variables.
#include "vtkCubeSource.h" #include "vtkPolyDataMapper.h" #include "vtkRenderer.h" #include "vtkRenderWindow.h" #include "vtkRenderWindowInteractor.h" #include "vtkInteractorStyleTrackballCamera.h" #include "vtkPolyData.h" #include "vtkGlyph3DMapper.h" vtkSmartPointer<vtkRenderer> renderer; vtkSmartPointer<vtkRenderWindow> renderWindow; vtkSmartPointer<vtkRenderWindowInteractor> renderWindowInteractor; vtkSmartPointer<vtkPoints> points;
Now we can write the function that will initialize the renderer and rendering window together with its interactor object.
void initializeRendering() {
renderer = vtkSmartPointer<vtkRenderer>::New(); //creating new renderer
renderer->SetBackground(1., 1., 1.); //setting background color to white
renderWindow = vtkSmartPointer<vtkRenderWindow>::New(); //creating renderer window
renderWindow->SetSize(1000, 1000); //setting size of the window
renderWindow->AddRenderer(renderer); //adding previously created renderer
renderWindowInteractor = vtkSmartPointer<vtkRenderWindowInteractor>::New(); //creating render
renderWindowInteractor->SetRenderWindow(renderWindow); //adding previously created renderer window
points = vtkSmartPointer<vtkPoints>::New(); //vector with points
}
I assume that we have our cell objects are stored in STL vector structure and that each cell contains information about its (x,y,z) position in 3D space. Having that we can write a really short function that adds our cells to the scene.
void addCells(vector<Cell>::iterator begin, vector<Cell>::iterator end) {
vector<Cell>::iterator cell = begin;
for (cell; cell<end; cell++)
points->InsertNextPoint(cell->x,cell->y,cell->z);
}
Finally, we can write the rendering function that is invoked in the main part of the code.
void visualizeLesion(vector<Cell>::iterator begin, vector<Cell>::iterator end) {
initializeRendering();
addCells(begin, end);
// Combine into a polydata
vtkSmartPointer<vtkPolyData> polydata = vtkSmartPointer<vtkPolyData>::New();
polydata->SetPoints(points);
vtkSmartPointer<vtkCubeSource> cubeSource = vtkSmartPointer<vtkCubeSource>::New();//you can use sphere if you want
vtkSmartPointer<vtkGlyph3D> glyph3D = vtkSmartPointer<vtkGlyph3D>::New();
glyph3D->SetSourceConnection(cubeSource->GetOutputPort());
glyph3D->ScalingOff();
glyph3D->SetInputData(polydata);
glyph3D->Update();
// Create a mapper and actor
vtkSmartPointer<vtkPolyDataMapper> mapper = vtkSmartPointer<vtkPolyDataMapper>::New();
mapper->SetInputConnection(glyph3D->GetOutputPort());
vtkSmartPointer<vtkActor> actor = vtkSmartPointer<vtkActor>::New();
actor->SetMapper(mapper);
// Add the actor to the scene
renderer->AddActor(actor);
// Render
renderWindow->Render();
vtkSmartPointer<vtkInteractorStyleTrackballCamera> style = vtkSmartPointer<vtkInteractorStyleTrackballCamera>::New(); //like paraview
renderWindowInteractor->SetInteractorStyle( style );
renderWindowInteractor->Start();
}
Snapshot of resulting rendering window is presented in Fig. 3 – and yes, I agree, it’s not very impressive. However, we can easily tweak it and make it significantly better.
Figure 3. Exemplary 3D tumor visualized using 3D glyph mapper in VTK.
First of all to have better a better perspective we can slice the tumor in half and take resulting parts slightly apart. This can be achieved by modification of the addCells function (note that in my setting center of the tumor is at (x,y,z) = (250,250,250)).
void addCells(vector<Cell>::iterator begin, vector<Cell>::iterator end) {
vector<Cell>::iterator cell = begin;
for (cell; cell<end; cell++) { if (cell->z < 250) { points->InsertNextPoint(cell->x,cell->y,cell->z);
} else {
points->InsertNextPoint(cell->x,cell->y,cell->z+50);
}
}
}
Fig. 4 shows the snapshot of the rendering window after slicing the tumor – much better, but there is still place for easy improvement. Let’s add some color!
Figure 4. Exemplary 3D tumor visualized using 3D glyph mapper in VTK.
I assume that each cell has also information about the number of free spots in its direct neighborhood. We will color the tumor based on that information.
First, we need to add additional headers and variables.
#include "vtkCubeSource.h" #include "vtkPolyDataMapper.h" #include "vtkRenderer.h" #include "vtkRenderWindow.h" #include "vtkRenderWindowInteractor.h" #include "vtkInteractorStyleTrackballCamera.h" #include "vtkPolyData.h" #include "vtkGlyph3DMapper.h" #include "vtkPointData.h" #include "vtkFloatArray.h" #include "vtkLookupTable.h" vtkSmartPointer<vtkRenderer> renderer; vtkSmartPointer<vtkRenderWindow> renderWindow; vtkSmartPointer<vtkRenderWindowInteractor> renderWindowInteractor; vtkSmartPointer<vtkPoints> points; vtkSmartPointer<vtkFloatArray> scalar;
Then we need to modify the function that initializes rendering
void initializeRendering() {
renderer = vtkSmartPointer<vtkRenderer>::New(); //creating new renderer
renderer->SetBackground(1., 1., 1.); //setting background color to white
renderWindow = vtkSmartPointer<vtkRenderWindow>::New(); //creating renderer window
renderWindow->SetSize(1000, 1000); //setting size of the window
renderWindow->AddRenderer(renderer); //adding previously created renderer
renderWindowInteractor = vtkSmartPointer<vtkRenderWindowInteractor>::New(); //creating render
renderWindowInteractor->SetRenderWindow(renderWindow); //adding previously created renderer window
points = vtkSmartPointer<vtkPoints>::New(); //vector with points
scalar = vtkSmartPointer<vtkFloatArray>::New(); //vector with scalar values of the points
scalar->SetNumberOfComponents(1); //one value per point
}
and addCells function
void addCells(vector<Cell>::iterator begin, vector<Cell>::iterator end) {
vector<Cell>::iterator cell = begin;
for (cell; cell<end; cell++) { if (cell->z < 250) { points->InsertNextPoint(cell->x,cell->y,cell->z);
} else {
points->InsertNextPoint(cell->x,cell->y,cell->z+50);
}
scalar->InsertNextValue((double)cell->nSpots/26.);
}
}
Finally, we can modify the visualization function.
void visualizeLesion(vector<Cell>::iterator begin, vector<Cell>::iterator end) {
initializeRendering();
addCells(begin, end);
//creating new color table
vtkSmartPointer<vtkLookupTable> colorLookupTable = vtkSmartPointer<vtkLookupTable>::New();
colorLookupTable->SetNumberOfColors(256);
colorLookupTable->SetHueRange( 0.667, 0.0);
colorLookupTable->Build();
// Combine into a polydata
vtkSmartPointer<vtkPolyData> polydata = vtkSmartPointer<vtkPolyData>::New();
polydata->SetPoints(points);
polydata->GetPointData()->SetScalars(scalar);
vtkSmartPointer<vtkCubeSource> cubeSource = vtkSmartPointer<vtkCubeSource>::New();//you can use sphere if you want
vtkSmartPointer<vtkGlyph3D> glyph3D = vtkSmartPointer<vtkGlyph3D>::New();
glyph3D->SetColorModeToColorByScalar();
glyph3D->SetSourceConnection(cubeSource->GetOutputPort());
glyph3D->ScalingOff();
glyph3D->SetInputData(polydata);
glyph3D->Update();
// Create a mapper and actor
vtkSmartPointer<vtkPolyDataMapper> mapper = vtkSmartPointer<vtkPolyDataMapper>::New();
mapper->SetInputConnection(glyph3D->GetOutputPort());
mapper->SetLookupTable(colorLookupTable);
vtkSmartPointer<vtkActor> actor = vtkSmartPointer<vtkActor>::New();
actor->SetMapper(mapper);
// Add the actor to the scene
renderer->AddActor(actor);
// Render
renderWindow->Render();
vtkSmartPointer<vtkInteractorStyleTrackballCamera> style = vtkSmartPointer<vtkInteractorStyleTrackballCamera>::New(); //like paraview
renderWindowInteractor->SetInteractorStyle( style );
renderWindowInteractor->Start();
}
Snapshot of sliced and colored tumor is presented in Fig. 5 – fast, good-looking and easy.
Figure 5. Exemplary 3D tumor visualized using 3D glyph mapper in VTK. Color represents the number of free spots in the direct neighborhood of each cell (jet colormap).
COMPUTE CANCER 
