Forking network generator is incorrect

[jmpf]

The forking network generator is incorrect. See 189d810 for the erroneous line.

Problem

This relates to Bernabeu et al. 2020. As described on page 3 of the supplementary information https://www.pnas.org/doi/suppl/10.1073/pnas.2007770117/suppl_file/pnas.2007770117.sapp.pdf Murray's law says that the vessel diameters scale with cube root. With fixed lambda (length/diameter) this means that the lengths also scale with cube root. In the cube the y-height of each bifurcation is taken to scale by 1/2 and the horizontal position is calculated from knowing the vessel length and vessel height. It's this horizontal position which is wrong. The actual vessel length in our forking networks does not match the description in the SI of the paper.

• vessel lengths are not as advertised so the geometry of all forking networks is wrong
• because vessel diameters are scaled correctly this means that lambda is not fixed to expected value
• because the Pries with memory solver uses a hard-coded parent vessel length (which is the expected length not the actual length) the CFL recovery is wrong
• All panels in Figure 4 of the paper have the potential to be wrong https://www.pnas.org/doi/10.1073/pnas.2007770117#fig04

Todo

• Assess how correcting the parent vessel length (to the actual length in these networks) will affect the Pries with memory solver in these incorrect geometries. See Fig 4b
• Assess how correcting the geometry will affect Fig 4a and Fig 4b. (How big are these length differences and are they noticeable on Fig4a?)
• Assess how correcting the geometry will affect Fig 4c and the SI data
• Potentially issue an erratum

[jmpf]

This what our current code produces for Figure 4b (before clipping and colour adjustments)

[jmpf]

See commit f74f08c
This Figure 4b with the same geometry but with the haematocrit solver using the actual vessel lengths (before clipping and colour adjustments)

There is a visible difference. In the original figure the oxygen concentration differences in the central region are more exaggerated. Perhaps this is because the expected vessel lengths are shorter than the actual length (meaning that the CFL recovers more in the new calculation).

[jmpf]

7135cec changes the default behaviour. The haematocrit converges to a slightly different place but, as can be seen here

and in ffea809, it's a negligible change to the original plot.

[jmpf]

See 1392ef5.

The correct geometry is now implemented and tested. This what the new code produces for Figure 4b (before clipping and colour adjustments).

It is noticeably less elongated than the original but has a really similar heterogeneous distribution is the oxygen field (see top of issue).

Next job is to clip out and colour adjust the central part of the figure so that it can be compared to Figure 4b. The squatness of the new geometry should be less noticeable when focusing solely on the central region of the network. After that, we can look at the violin plots.

[jmpf]

Summary

Figure 4b difference

This is my best shot at reproducing the heat map tissue oxygenation. The axes labels and colours are slightly off in both the "original" and "corrected". The "corrected" version is not so wide but the colour trends are practically the same.

Original	Corrected

Figure 4c difference

The violin plots are different if you flick between them or if you compare pixels (below). They are qualitatively the same.

Original	Corrected

Pixel-by-pixel comparison

Original Figure 4 from PNAS paper (for comparison)

Notes

Figure 4a correction is like Fig. 4b but with uniform oxygen.

Similar corrections could be done to the SI. However these corrections can be imagined because it is noted that:

• Figure S5 shows how Fig. 4a and 4b fit into a larger network (so the correction will be a simple horizontal re-scaling);
• Figure S6 shows the same data as the violin plots but with box/whiskers, means and standard deviations.