Surface-guided computing to analyze subcellular morphology and membrane-associated signals in 3D

Signal transduction and cell function are governed by the spatiotemporal organization of membrane-associated molecules. Despite significant advances in visualizing molecular distributions by 3D light microscopy, cell biologists still have limited quantitative understanding of the processes implicated in the regulation of molecular signals at the whole cell scale. In particular, complex and transient cell surface morphologies challenge the complete sampling of cell geometry, membrane-associated molecular concentration and activity and the computing of meaningful parameters such as the cofluctuation between morphology and signals. Here, we introduce u-Unwrap3D, a framework to remap arbitrarily complex 3D cell surfaces and membrane-associated signals into equivalent lower dimensional representations. The mappings are bidirectional, allowing the application of image processing operations in the data representation best suited for the task and to subsequently present the results in any of the other representations, including the original 3D cell surface. Leveraging this surface-guided computing paradigm, we track segmented surface motifs in 2D to quantify the recruitment of Septin polymers by blebbing events; we quantify actin enrichment in peripheral ruffles; and we measure the speed of ruffle movement along topographically complex cell surfaces. Thus, u-Unwrap3D provides access to spatiotemporal analyses of cell biological parameters on unconstrained 3D surface geometries and signals.


150
In summary, u-Unwrap3D generates bijective mappings of a given genus-X surface between 5 equivalent 151 surface representations; Cartesian 3D, ( , , ), genus-0 reference 3D, ref ( , , ), the unit 3D sphere, 152 2 ( , , ) , topographic 3D, ( , , ) , and the 2D plane image, ( , ) (Fig. 2a), while simultaneously 153 transforming Cartesian 3D to topographic 3D volumes (Fig. 2b). This was made possible by two crucial 154 choices; the use of cMCF and voxelization to construct a genus-0 reference surface, ref ( , , ) to realise 155 spherical parameterization (Step 1) and the implementation of an efficient numerical scheme to relax area 156 distortion on the 3D sphere (Step 3). The former allows us to construct ref ( , , ) as a proxy of the genus-157 X ( , , ) surface mesh and to unwrap this 3D surface into one 2D ( , ) image, instead of requiring multiple 158 2D ( , ) images, which simplifies downstream analysis [56][57][58][59] . The latter ensures that the unwrapped 2D ( , ) 159 image captures the salient surface features of ref ( , , ), and by extension the genus-X ( , , ) surface. 160 Importantly, the bijectivity of the mappings guarantees that for any point on any of the surface or volume 161 representations matching points exist on any of the other surfaces or volumes. Moreover, the bijectivity 162 guarantees preservation of the point topology, i.e. a series of points ordered in clockwise fashion on one 163 surface representation maps to a series of points ordered in the same way on any of the other surface 164 representations and preserves the local neighbourhood relationships. As a result, we can apply mathematical 165 operations defined in any one of the representations and map the results to any other. u-Unwrap3D thus 166 supports the optimal spatiotemporal analysis of unconstrained surface geometries and associated signals.

168
Validation of u-Unwrap3D on diverse surface motifs 169 We validated the generality and performance of u-Unwrap3D by application to 66 single cell images acquired 170 by high-resolution light sheet imaging 45,60 . The dataset span morphologically diverse cells with blebs, 171 lamellipodia and filopodia. The cell surfaces were meshed with marching cubes and segmented within the u-172 Shape3D software 45 . Small errors in the initial segmentation and meshing process cause high-order genus 173 surfaces with topological holes and handles, which cannot be unwrapped directly (Extended Fig. 2a). Holes 174 can also generate non-watertight surface meshessurfaces that are not closed and have no clearly defined 175 inside volume 48,61 possessing potentially complex internal volumetric structures that violate the assumptions 176 of standard 3D mesh processing algorithms. 177 178 We first tested the number of input cell surfaces for which u-Unwrap3D could successfully run all steps 1-6 179 and compute all 5 of the representations as a measure of generality and robustness. Notably only 6/66 (11%) 180 input cell surfaces were genus-0 and only 36/66 (55%) were watertight (Extended Fig 2b). In 63/66 cases, 181 (>95%) we successfully ran all steps and obtained all representations (Extended Fig 2b). The three failures 182 occurred in scenarios, in which the holes and handles remaining after the application of cMCF were still too 183 large for the volume dilation to generate a genus-0 reference surface after remeshing (Fig. 1b, Step 1) 184 (Extended Fig. 2c). In all successful cases, cMCF and binary voxelization under volume dilation generated 185 genus-0 reference surfaces within a median of 10 iterations (Extended Fig 2c, c.f. lamellipodia). Fig. 2c  186 shows extracted representations for challenging examples with blebs, lamellipodia and filopodia (Suppl. 187 Video 3-5). 188 189 We next tested the robustness and performance of the ref ( , , ) spherical parameterizations, (Fig. 1b,  190 Steps 2-3). Extended Fig. 3a confirms that the quasi-conformal spherical parameterization (Step 2) minimizes 191 the conformal error to the ideal value of 1, with the largest error in cells with filopodia (1.016±0.013). We also 192 verified the need to relax local area distortion. Whilst quasi-conformal spherical parameterization Ω 2 ( , , ) 193 is equiareal for blebs (0.978 ± 0.037), the median area distortion showed that the surface fraction of 194 lamellipodia was down to 0.432 and in filopodia to just 0.140 with respect to their original area fraction on the 195 reference 3D, ref ( , , ) surface, let alone ( , , ). Our scheme for area distortion relaxation (Methods) 196 produces a quasi-equiareal spherical parameterization (Step 3) in blebs (0.985±0.012), and successfully 197 achieves the ideal value of 1 in lamellipodia (1.000±0.000) and filopodia(1.000±0.001), within a maximum 198 median of 23 iterations for lamellipodia (Extended Fig.3b, Table i). We further tested the ability of our 199 relaxation scheme to balance the trade-off between the two extremes of conformal to equiareal spherical 200 parameterizations using different stopping criteria (Extended Fig.3b, Table ii-iv). The initial parameterization 201 without any area-distortion relaxation (iteration 0) is by design conformal but also found to satisfy the most 202 isometric parameterization (MIP) 62 . Running for Ω < a maximum of 50 iterations yields an equiareal 203 parameterization for all motifs. At ≈ 1 2 Ω iterations the relaxed mesh jointly minimizes the summation ( + 204 ln , Methods) of conformal ( ) and area distortion ( ) errors. At ≲ Ω iterations the relaxed mesh is the 205 area-preserving MIP 63 . As expected, this latter parameterization does not fully minimize area distortion (blebs 206 (0.997±0.011), lamellipodia (0.979±0.005) and filopodia (0.980±0.028)) but exhibits slightly lower conformal 207 errors and consequently higher quality faces than a pure equiareal mapping. 208 209 Lastly, we tested how accurately topo ( , , ), which defines the topographic 3D mesh, ( , , ) (step 5) 210 mapped back into Cartesian coordinates reconstructs the input surface, ( , , ). For all cells, ( , , ) was 211 computed with a ( , ) image grid size of 1024x512 pixels. The aspect ratio, 2 x ( = 512) was chosen 212 to preserve the ratio between the equatorial circumference and the length of the arc between north and south 213 poles of a sphere. Compared to the input surface ( , , ), ( , , ) is lower genus and provides higher face 214 quality (Extended Fig. 3c). We assessed the discrepancy between topo ( , , ) and ( , , ) using 4 metrics; 215 Chamfer distance (CD), sliced Wasserstein distance ( 1 ) 64 , and differences in total surface area (Δ ) and 216 volume (Δ ) (Extended Fig. 3d, Methods). Considering inevitable rasterization errors when mapping the 217 floating-point precision 3D sphere Ω 2 ( , , ) to ( , ) defined on an integer ( , ) image grid, we measured 218 low vertex position errors according to CD and 1 . Cells with lamellipodia had the lowest error (median 219 CD=1.77 voxel, 1 =0.93 voxel) and cells with blebs were slightly worse (median CD=2.79 voxel, 1 =4.28 220 voxel), likely due to their intrinsically small height (small topographic ). As one would expect, cells with long, 221 thin filopodia displayed the largest discrepancies (median CD=10.33 voxel, 1 =18.45 voxel). 222 Correspondingly we measured a small Δ (+1.2%) and Δ (+7.9%) for cells with lamellipodia. Δ was larger 223 for cells with blebs (-11.6%) and measured to be too large for filopodia (-55.3%) when compared to Δ 224 differences measured after making ( , , ) watertight (+4.2% blebs, +3.4% filopodia). Visualization of 225 exemplar cells show good geometric correspondence between topo ( , , ) and ( , , ) (Extended Fig. 3e).

226
Salient surface features were largely captured, albeit smoothened and blurred in topo ( , , ) when local 227 surface regions were underrepresented due to being distant relative to ref ( , , ) (Extended Fig. 3e, black 228 triangles, 1 st row blebs and 4 th row lamellipodia). Most of the primary morphological features, namely the 229 length and thickness of long, thin filopodia (Extended Fig. 3e, green triangles), except those located both 230 densely together and distant relative to ref ( , , ) (Extended Fig. 3e, red triangles), were captured. u-231 Unwrap3D was able to capture both the complex lamellipodia folds and curved cell bodies to high accuracy 232 (Extended Fig. 3e, row 2,4). Closer inspection of topo ( , , ) and ( , , ) in these cells traced a large Δ 233 to meshing errors in the input surface ( , , ), which caused internal volumetric structures to be merged 234 into the cell surface representation, and overestimation of total surface area. These errors affect the CD and 235 1 to lesser extent.

237
In summary, our results demonstrate that u-Unwrap3D is robust and applicable to process unconstrained 238 geometries. For maximum resolution of high curvature surface features, a genus-0 reference surface 239 ref ( , , ) proximal to the input surface ( , , ) is recommended with a large ( , ) grid size and small 240 step sizes when propagating ref ( , ). However, these choices depend on the quality of the input surface 241 mesh ( , , ), which depends on the robustness of cell segmentation in the face of noisy image raw data. 242 243 u-Unwrap3D enables unsupervised instance segmentation of subcellular surface motifs 244 The We demonstrate the segmentation of individual protrusion instances, capturing motifs identified by uShape3D, 282 but without the need for training annotations. By construction, in topographic 3D ( , , ) space all surface 283 protrusions are oriented upwards with increasing . Moreover, the tops of protrusions are individually 284 separated as local regions of high topographic mean curvature, which we identify by thresholding and 285 applying connected component labelling in the topographic volume, ( , , ) (Fig. 3c). Mapping the 286 segmented regions back onto the topographic mesh ( , , ), we diffuse these initial 'seed' labels across 287 the surface using a combined geodesic distance and dihedral angle affinity matrix to naturally segment the 288 'stem' of the individual protrusions (Methods). The dihedral angle measures the discontinuity in local mean 289 curvature. It incorporates the prior intuition that the boundaries of a label should expand faster on local 290 surfaces of homogeneous curvature such as that on a 'hill', compared to another label experiencing large 291 curvature differences in its local surface region, such as in a valley between multiple 'hills'. The combined 292 affinity matrix thus introduces a morphology-aware competition between segmentation labels and provides a 293 biophysical rationale for defining which surface patches belong to individual 'seed' protrusions. Furthermore, 294 the dihedral angle is large between a protrusion and the main cortical cell surface and thus serves as a soft 295 stopping criteria for diffusion (Extended Fig. 4d) in addition to applying the binary protrusion segmentation 296 from above. Lastly, we filter out protrusions that are too small and close any small holes using the Cartesian 297 3D surface area. The final segmentation result qualitatively and quantitatively agrees with that obtained by 298 supervised u-Shape3D for lamellipodia ( Fig. 3d) but yields more contiguous labels and is less prone to over-299 segmentation (Fig. 3e). Importantly, this segmentation strategy is applicable, even when not all protrusions 300 are equiareally represented in ( , , ), as shown by the segmentation of the majority of blebs and filopodia 301 in exemplar cells (Extended Fig. 4d).

303
Individual segmented protrusions are genus-0 open-surface 3D submeshes that can be directly mapped to 304 the 2D plane (Fig. 3f). This allows us to further refine the segmentation, for example, by detecting and splitting 305 under-segmented blebs by a gradient watershed algorithm (Methods). Thanks to bijectivity, the refined 306 segmentation labels can be transferred back onto the original surface mesh (Fig. 3f,g, c.f. before and after 307 refine, black triangles). Both the before (adjusted normalized mutual information, NMI=0.57) and the after 308 refinement (adjusted NMI=0.54) segmentations agree with u-Shape3D. However, u-Unwrap3D segmented 309 blebs are more complete, with more blebs of larger surface area (150 blebs in total). In contrast, u-Shape3D 310 over-segments small blebs (742 blebs) that were found to originate from erroneous meshing of internal 311 structures in ( , , ) (see Extended Fig.3e). This example illustrates the potential pitfalls of identifying motifs 312 from local surface patches only with potentially imperfect input 3D meshes. In contrast, with u-Unwrap3D any 313 . CC-BY 4.0 International license made available under a (which was not certified by peer review) is the author/funder, who has granted bioRxiv a license to display the preprint in perpetuity. It is The copyright holder for this preprint this version posted April 20, 2023. ; surfaces internal to the cell volume are readily removed when mapped into topography as these surfaces 314 have -coordinates less than the reference surface, ref ( ref = smooth ( , ), , ).

316
Finally, the representation ( , , ) enables partitioning of the input cell volume into the sum of a reference 317 volume representing the underlying cortical cell body and the unique volume occupied by individual 318 protrusions (Fig. 3h). We do so by ( , )-parameterizing the reference cortical surface submesh, ref ( , , )  319 after removing all individual protrusion submeshes as a grayscale image such that the pixel intensity value 320 at ( , ) is the respective -coordinate, and inpainting the missing -coordinates at ( , ) coordinates 321 corresponding to the subtracted protrusions to generate a full reference binary volume (Extended Fig. 4e, 322 Methods). For the protrusions, we first devised a marker-controlled lateral watershed depth propagation to 323 diffuse the surface-based protrusion segmentation labels uniquely throughout the full topographic ( , , ) 324 space slice-by-slice, top-to-bottom (Extended Fig. 4f). Then, individual protrusion volumes were generated 325 by masking the propagated label volume with the reference binary volume (Extended Fig. 4g). We compared 326 our topography-guided decomposition strategy to a fully Cartesian 3D mesh processing approach whereby 327 individual protrusion submeshes were first closed by constructing a surface patch that minimized the local 328 bending (or harmonic) energy 72 (Methods). The closed reference volume was then generated using all such 329 patches to impute the residual holes in the raw reference surface. Fig. 3i  However, the 3D mesh processing protrusions consistently under-measure larger protrusions. Crucially, the 332 imputed reference volume appears artefactual. Where protrusions were located, the surface is overly smooth, 333 and even involuted. These regions contrast starkly with non-imputed surface areas between protrusions, 334 creating artificial 'peaks' and 'ridges' of high mean curvature (black arrows). In comparison, the topography-335 guided reference volume is mechanically more plausible.

337
In summary, the ability to map freely between topographic and Cartesian 3D surfaces and their respective 338 volumetric representations enabled us to design simple and generalizable methods to detect and segment in 339 an unsupervised fashion individual morphological motifs from unconstrained surface geometries. To 340 bijectively map the topographically segmented surface protrusion labels onto the 2D plane, which is an 341 optimal representation for tracking the segmented motifs, we developed a topographic cMCF for u-Unwrap3D 342 (Fig. 3j

Individual blebs dynamically recruit Septins to local surface regions during retraction 357
We first analysed potential relations between dynamic surface blebbing and the recruitment of Septins. Blebs 358 are globular membrane protrusions of 1-2 μm diameter that are thought to extend in areas of localized 359 membrane detachment from the actin cortex 75,76 . Intracellular pressure expands the budding blebs outward, 360 followed by rapid assembly of a contractile actomyosin network that yields retraction. Cycles of protrusion 361 and retraction have been described to last a few tens of seconds. Associated with the cycles are molecular 362 activities both driving and responding to the morphological dynamics. One such process is the assembly of 363 Septin protein polymers at sites of negative curvature (from a cell-external perspective) emerging at the bleb 364 necks. Our previous work 38 has shown that disrupting the bleb cycle diminishes Septin assembly at the cell 365 surface. Here we now exploit the ability of u-Unwrap3D to track individual bleb cycles and quantify Septin 366 accumulation by remapping surface morphology and a fluorescent marker of Septins to an appropriate 2D 367 representation. We acquired 3D volumes of SEPT6-GFP-expressing MV3 melanoma cells every 1.2s for 200 368 timepoints. As the cortical cell body exhibits little temporal variation and blebs protrude normally to the surface, 369 the temporal mean cell surface, ̅ ( , , ) is a good proxy of the cell cortex. We apply u-Unwrap3D to ̅ ( , , ) 370 to create a common static ( , , ) coordinate space for computing topographic 3D ( , , ) and 2D planar 371 ( , ) representations for each timepoint. ̅ ( , , ) was computed by meshing the mean binary volume 372 across individual binary voxelizations of the cell surface over all 200 timepoints (Methods). Note the 373 construction of a common topographic space from a single mesh for a timelapse is computationally efficient 374 but applicable only if cell shape changes lie within the Cartesian subvolume mapped by the topographic 375 space. For large shape changes u-Unwrap3D should be applied to individual timepoints and spatiotemporal 376 registration used to align surfaces to a common reference using ref ( , , ) or ( , ) representations. The 377 segmentation tools discussed above were used to detect all bleb instances from ( , , ) and mapped to 378 ( , , ) and ( , ) for each timepoint (Fig. 4a, Suppl. Video 7). Similarly, the computed mean curvature 379 and the normalized Septin intensity surface signals ( ( , , )) from Cartesian 3D were mapped to 380 topographic 3D ( ( , , )) and into 2D, ( ( , )), to enable simple bleb tracking and timeseries analysis. 381

382
To track blebs in 2D we computed the bounding box of individual bleb instances in every timepoint after 383 appropriate image padding to account for spherical periodicity (Fig. 4b left). In ( , ) bleb dynamics can 384 readily be followed by an established 2D multi-object bounding box tracker 77 with mean curvature optical 385 flow-guided bipartite matching ( Fig. 4b  portion of the cell can be maintained in-focus. Again, u-Unwrap3D's bijectivity between 3D and 2D 389 representation enabled us to easily map a manually annotated out-of-focus subvolume onto ( , , ) and 390 into ( , ) to retain for analysis only the bleb tracks that remain within the in-focus surface regions (Fig. 4c).

391
The distribution of individual bleb track lifetimes showed a peak at 14s and a long tail up to 240s, suggestive 392 of a mixture of short-and long-lived blebs (Fig. 4d). The mean bleb lifetime of 27s corresponded well with a 393 30s periodicity given by the first peak of the temporal autocorrelation of mean curvature ( ( , , , )) in 394 Cartesian 3D (Fig. 4e). This validates at the population level the accuracy of single bleb tracking after 395 projection and segmentation of mean surface curvature in 2D. 396 The temporal autocorrelation curves of mean curvature and Septin intensity showed a high level of similarity, 397 suggesting co-fluctuation of the two surface signals. Indeed, we had previously shown that surface regions 398 of high Septin intensity with negative surface curvature for at least 30s display a correlation between negative 399 curvature value and Septin intensity 38 . Whilst the majority of Septin pulses endured only one cycle of bleb 400 formation and retraction, de novo formation of stable Septin structures appeared to be driven by several 401 Septin pulses occurring in short succession. We thus hypothesized that these were formed by iterative bleb-402 driven curvature generation events resulting in local levels of Septin oligomers surpassing a threshold 403 necessary for inter-oligomer polymerization and enabling stabilization through formation of higher-order 404 structures 38 . The ability to spatiotemporally track individual blebs enabled us now to quantitatively test this 405 model. For the duration of each tracked bleb, we sampled within the 2D bounding box distortion-corrected 406 timeseries of bleb surface area, on-/off-bleb surface mean curvature and Septin intensity (Methods). We 407 used the 2D ( , ) bounding box to define a bleb's spatial area-of-influence and its surface area. The 408 Cartesian 3D ( , , ) box area was taken as the bleb's 3D surface area in each tracked frame. Within the 409 2D bleb bounding box, 'on-bleb' is the largest spatially contiguous region of high mean curvature. The 410 remainder area is 'off-bleb'. A single curvature threshold was computed by 3-class Otsu thresholding over all 411 ( ( , , )) to define the regions of low/high mean curvature in ( , , ). Using the extracted timeseries we 412 reconstructed the temporal profile of bleb area, on-/off-mean curvature and Septin intensity of a single mean 413 bleb event in a window of 35s centered on the timepoint of maximum bleb area averaged over 545 single 414 bleb events from 480 bleb tracks. We sample ±17.5s before and after the timepoint of maximum bleb area 415 to exceed the 30s periodicity inferred from temporal autocorrelation by 5s and capture the full dynamics. We 416 then broke the mean timeseries into four distinct temporal phases of bleb-mediated Septin recruitment, each 417 ≈5s long (Fig.4f labels A-D). In phase A, the bleb begins to expand, increasing its surface area, accompanied 418 by a sharp increase in on-bleb and a decrease in Septin intensity as the plasma membrane detaches from 419 the actin cortex. The expansion also reduces off-bleb , causing a decrease in Septin intensity off-bleb, 420 presumably due to disrupting Septin structures. In phase B, the bleb reaches maximum size and then begins 421 to retract with decreasing surface area. Interestingly, unlike mean curvature, which begins to decrease before 422 the maximum bleb area, the change in area is symmetrical, occurring ±2.5s relative to the time of maximum 423 bleb area. Coincident with the bleb increasing to a maximum area, off-bleb decreases to a minimum and 424 Septin intensity both on/off bleb stabilizes at a minimum. As the bleb retracts and the actin cytoskeleton 425 reassembles, off-bleb increases and Septin intensity begins to increase both on/off bleb. In phase C (+2.5s 426 to +8.0s after peak bleb area), the bleb area continues to decrease but at a slower rate. Unexpectedly, the 427 off-bleb negative curvature plateaus at a value lower than its starting value (before phase A) instead of 428 continuing to increase. Concurrently, Septin intensity undergoes the greatest rate of increase on both on-429 bleb and off-bleb surfaces such that at the end of phase C, the Septin intensity is at the levels before phase 430 A (c.f. -15s to -10s). In phase D (+8.0s to +14.0s), bleb area and on/off-bleb recover to pre-expansion 431 levels while Septin intensity continues to increase before plateauing on both on/off-bleb surfaces. Beyond 432 phase D, the next bleb cycle begins, with similar temporal characteristics to phase A. Altogether these results 433 indicate that blebs generate optimal curvature dynamics during retraction to recruit Septin polymers to the 434 surface regions around blebs. Moreover, the data support our model of Septins being recruited to negative 435 curvature patches in a cyclic fashion, where each bleb formation-retraction drives the local accumulation of 436 a few more oligomers until a threshold concentration is reached to trigger inter-oligomer polymerization of a 437 stable Septin assembly. Notably, none of these observations could have been made without u-Unwrap3D. cover glass tilted at ≈45 0 (Extended Fig. 5c, left). Direct 3D plane-fitting to ( , , ) to determine a precise 478 angle is sensitive to outlier points and shapes deviating from an elongated ellipsoid. Thanks to the ( , ) 479 . CC-BY 4.0 International license made available under a (which was not certified by peer review) is the author/funder, who has granted bioRxiv a license to display the preprint in perpetuity. It is The copyright holder for this preprint this version posted April 20, 2023. ; https://doi.org/10.1101/2023.04.12.536640 doi: bioRxiv preprint representation of u-Unwrap3D mapping the curved proximal cell surface and the planar cell bottom to the 480 upper and lower half of the unwrapped 2D image respectively, we could readily annotate in 2D ( , ) the cell 481 bottom and fit a 3D plane to only this surface patch in ( , , ) (Extended Fig. 5c, middle, right). Similarly, we 482 could gate ROI tracks in ( , ) and compute the lateral speed only for those associated the lamellipodia and 483 lamella surface. Doing so we found two populations in the speed histogram (Fig. 5d). Visualizing the speed 484 on ( , , ), the faster population corresponds to higher curvature ruffles with speeds ranging from 2-10 485 μm/min, and an average speed of 4.2 μm/min (Suppl. Video 8). This is at least two times faster than the 486 slower population surface retrograde actin flow ranging from 0-1 μm/min, which are consistent with the flow 487 speeds we used to measure by 2D fluorescent speckle microscopy in the lamella of epithelial cells 85  Analyzing the spatiotemporal organization of molecular distributions and signaling activities on cell surfaces 503 in 3D has been limited by the lack of methods to represent, track and process these dynamics. Here we 504 introduce a surface-guided computing framework, referred to as u-Unwrap3D, to bijectively map a genus-X 505 Cartesian 3D surface to equivalent surface and volume representations that are optimally suited for a distinct 506 analytical task. The mappings rely on two critical insights: i) the engineered surface deformation of ( , , ) 507 to generate a genus-0 ref ( , , ) for which a 3D spherical parameterization exists; ii) a novel, efficient 508 algorithm to relax geometric distortion on the 3D sphere in a bijective and tunable manner. Insight i) is applied (>90% of cases) in a manner robust to the input surface mesh quality and that it accurately captures 513 the cell geometry to transfer salient surface features, morphological or molecular, between all representations. 514 We note that this 90% is an underestimation of the applicability. The validation dataset was assembled to be 515 deliberately heterogeneous and was not segmented with the downstream aim of surface mapping. In practice, 516 segmentation algorithms continue to improve and postprocessing techniques can be used to create improved 517 surface meshes. Moreover we can leverage surface meshing algorithms that are more sophisticated than 518 marching cubes, such as dual contouring 86 and shrink-wrapping 87 , to guarantee watertight mesh creation. 519 For timelapse acquisitions, as we showed in Fig. 4,5, the situation is even simpler. Only one timepoint or 520 average surface is needed to generate a single common ̅ ref ( , , ) to unwrap all timepoints. u-Unwrap3D 521 puts in place a generic platform for the spatiotemporal processing of unconstrained cell geometries. 522 u-Unwrap3D is applicable for arbitrary genus-X 3D surfaces, ( , , ) wherever a genus-0 reference surface, Whilst we can increase the range of morphological hole closing in the latter case, we restrict dilation to a 528 maximum 3-5 voxels to minimise smoothing out protrusion features in ref ( , , ). Mesh surgery methods 529 have been developed in computer graphics to make non-watertight 3D meshes watertight 88,89 , and fix 530 imperfections like holes and handles to reduce genus 90,91 to generate higher quality 3D meshes. These have 531 yet to be applied and fully tested for complex cell surfaces. Future development will investigate how to apply 532 such procedures to allow u-Unwrap3D to be applied to input meshes of any quality. In our data we have paid 533 attention during the acquisition to generating sufficient foreground-to-background contrast for reliable surface 534 segmentation, thus minimizing mesh defects. 535 . CC-BY 4.0 International license made available under a (which was not certified by peer review) is the author/funder, who has granted bioRxiv a license to display the preprint in perpetuity. It is The copyright holder for this preprint this version posted April 20, 2023. ; https://doi.org/10.1101/2023.04.12.536640 doi: bioRxiv preprint The concept of geometrical reduction of 3D into 2D geometry through the choice of an optimal coordinate 536 transformation has long existed in mathematics and physics to simplify mathematical manipulation and 537 plotting. For example, parametric coordinates describe the sphere, the cylinder, Mobius strip and helicoids 538 amongst others 92 . In computer graphics this is realised in the common practice of mapping between surfaces 539 through simpler intermediaries: for example, the texture mapping of arbitrary surfaces by optimal surface 540 cutting and mapping of the cuts into individual 2D shapes 57 and ( , ) surface parameterization by cutting 541 and gluing individually mapped 2D planar patches 93 or mapping to canonical shapes such as the triangle 94 , 542 plane 39 or polyhedra 95 to minimize distortion. u-Unwrap3D draws inspiration from this thinking. Through the 543 availability and development of rationalized multiple 3D-to-2D representations, u-Unwrap3D projects 544 analyses that would otherwise require specialised mathematical operations into a sequence of simpler, 545 computationally tractable procedures for 3D mesh processing, image processing and machine learning with 546 specific consideration for single cell biology. Unlike computer graphics benchmarks, surface protrusions are 547 irregular, high curvature and dense. First, we specifically chose representations that map the whole cell 548 surface with well-behaved topologies such as the sphere and plane and designed a relaxation scheme to 549 guarantee interpolation between the minimal conformal and area distortion. This bypasses the numerical the interplay of dynamic membrane morphology and associated signals with more volumetric 556 nuclear/cytoplasmic signalling. Lastly, we designed mappings to underlying representations, including 557 Cartesian 3D, 3D sphere and 2D plane that are standard inputs in computer vision and machine learning. 558 With u-Unwrap3D standard computer vision, machine learning methods become directly applicable to 559 computing tasks on rugged complex surfaces. More recently, research into combining different geometric 560 representations are state-of-the-art in addressing complex computational problems such as multiple 2D 561 image views to inform 3D mesh reconstruction 96,97 , or 3D mesh vertex coordinates with 2D unwrapped images 562 for feature extraction 98,99 . u-Unwrap3D is fully complementary to these research developments -unifying the 563 different representations into a single surface-guided computing framework for downstream analysis. The u-564 Unwrap3D framework is made available as a Python library. The resources and validation provided by this 565 work will aid the cell biology community to generate testable hypotheses of the spatiotemporal organization 566 and regulation of subcellular geometry and molecular activity. 567

568
Acknowledgments 569 We thank Hanieh Mazloom Farsibaf and Zach Marin for discussions on mesh processing. We also thank 570 Qiongjing Zou for putting the Python library into a continuous integration coding framework. Funding for this 571 work in the Danuser lab was provided by the grant R35 GM136428 (NIH). AW is a fellow of the Jane Coffin 572 Childs Memorial Fund. GMG is a Damon Runyon Fellow (DRG 2422-21). Competing Interests 581 The authors declare no competing interests. 582 583 Data Availability 584 All data presented herein are available from the corresponding author upon request. 585 586 Code Availability 587 u-Unwrap3D is open source, developed as a Python library and is available at: 588 https://github.com/DanuserLab/u-unwrap3D. 589 . CC-BY 4.0 International license made available under a (which was not certified by peer review) is the author/funder, who has granted bioRxiv a license to display the preprint in perpetuity. It is The copyright holder for this preprint this version posted April 20, 2023. ; https://doi.org/10.1101/2023.04.12.536640 doi: bioRxiv preprint 590 591 Methods 592 u-Unwrap3D framework 593 Following we describe the algorithms underpinning each step of u-Unwrap3D depicted in Fig.1. 594 Step 1: ( , , ) to ref ( , , ) 595 Conformalized mean curvature flow (cMCF). cMCF 49 modifies the mean curvature flow (MCF) to avoid the 596 formation of pinches and collapsed vertices that compromise bijectivity and cause early flow termination for 597 watertight meshes with high curvature features. We find cMCF also reduces the size of holes and handles in 598 Beltrami is the divergence of the gradient, Δ = ∇ ⋅ ∇ with respect to the local mesh metric to get 607 are called the mass ( ) and Laplacian ( ) matrices at time, . 609 Substituting this notation and rearranging, the linear algebra MCF equation is 610 ( + ) is then computed from ( ) by direct matrix inversion. cMCF modifies the MCF equation above by 612 using the Laplacian matrix at time = 0, 0 for all timepoints. The Laplacian matrix is a measure of stiffness 613 between local mesh faces, see active contour cMCF below. Using 0 for all timepoints instead of recomputing 614 implicitly constrains mesh faces to retain the same aspect ratio and this conformalizes the flow. 615 We use the libigl library 101 with the cotangent Laplacian and barycentric mass matrix as the default 617 implementations of , respectively. We improve the numerics of solving cMCF by normalization of the 618 surface area and recentering of vertex coordinates at the origin after each iteration as recommended in Alec is the cortical cell shape without protrusions. We find that this corresponds to finding the 'elbow point' in the 625 mean absolute Gaussian curvature (Fig. 1a) and not the convergence limit of cMCF which is the sphere 49 . 626 The difference in the mean absolute Gaussian curvature over vertices, Δ| | ̅̅̅̅̅ = | | ̅̅̅̅̅ − | −1 | ̅̅̅̅̅̅̅̅ between 627 successive iterations − 1 and is used as an automatic stopping criterion for cMCF, stop = max ( min , ). 628 min is a user-specified minimum iteration number and is the first iteration for which Δ| | ̅̅̅̅̅ exceeds a user-629 specified threshold, Δ| | ̅̅̅̅̅ > Δ ℎ ℎ . We compute the discrete Gaussian curvature, at a vertex , given quasi-conformal). In practice, we found conformal errors = 0 (Extended Fig. 3). 654

742
Mesh displacement by active contour cMCF 743 Active contours, or 'snakes' 108 , define the boundary of an image region by minimizing its contour energy, .

744
The contour energy is the sum of an internal, and an external energy, , diffeomorphically. We refer to this as active contour cMCF in this paper. To evolve mesh vertices, ( ) normal 758 to the surface in equal -steps, we set = Φ, the signed distance function, and solve ( − ) ( + 1) = 759 ( ( ) + ∇Φ) iteratively, with ∇Φ evaluated at ( ) for each iteration. A positive moves a cell surface 760 mesh normally outwards from the cell and a negative moves the mesh normally into the cell.

762
Quantification of geometric deformation errors for meshes 763 Surface mappings do not conserve local geometrical measures like angles, edge lengths and face area. 764 Quantification of the distortion in these measures enables a task-specific optimization of the mapping and 765 correction of statistical measurements made on the mapped surface. There are two primary geometric 766 distortions to quantify; conformal and area distortion error (Extended Fig. 1). An isometric deformation is one 767 with no error; both conformal and area distortion errors are 0. 768 . CC-BY 4.0 International license made available under a (which was not certified by peer review) is the author/funder, who has granted bioRxiv a license to display the preprint in perpetuity. It is The copyright holder for this preprint this version posted April 20, 2023. ; https://doi.org/10.1101/2023.04.12.536640 doi: bioRxiv preprint Conformal error. The conformal or quasi-conformal error, measures the extent the shape of a mesh 770 element , e.g. a triangle face, is stretched. It is 0 if the relative distances between vertices and the angles 771 between edges are preserved after the mapping. We compute of mapping triangle Δ to Δ in 3D by 772 first isometrically projecting all triangles into 2D. Let = ( 1 , 1 , 1 ), = ( 2 , 2 , 2 ), = ( 3 , 3  Stopping criteria for area distortion relaxation of 2 ( , , ) 807 We used three additional stopping criteria to demonstrate intermediate area distortion relaxation between 808 fully conformal, 2 ( , , ) and fully equiareal, Ω 2 ( , , ) parameterization. We use the same nomenclature 809 as for the above discussed geometric deformation error for meshes.

811
Most isometric parametrization (MIPS) error. The MIPS 62 error is defined 2 1 + 1 2 and is minimal when 1 = 812 2 . This error is trivially minimal for a conformal spherical parametrization ( 1 2 = 1) (Extended Fig. 3b). 813 814 . CC-BY 4.0 International license made available under a (which was not certified by peer review) is the author/funder, who has granted bioRxiv a license to display the preprint in perpetuity. It is ) . We use this metric 815 with = 1, which measures the uniformity of stretch distortion over a surface. This error yields near-equiareal 816 spherical parametrization (Extended Fig. 3b). 817 818 Isometric error. We observed tradeoff of conformal error, and log area distortion, log on a similar 819 magnitude scale such that their summation has a unique global minima (Extended Fig. 1e). We thus define 820 an isometric error metric, (1 − ) 1 2 + ( ) log with a constant ∈ [0,1] to weight the relative importance of 821 jointly minimizing conformal and area distortion error. We use = 0.5 in Extended Fig. 3b.

823
Assessment of geometrical difference between two meshes 824 825 Four metrics were used to assess the difference between two meshes 1 and 2 possessing different number 826 of vertices and faces; chamfer distance (CD), Wasserstein-1 distance ( 1 ), the difference in surface area 827 (Δ ) and the difference in volume (Δ ). CD is the mean Euclidean distance between all vertices of 1 when 828 matched to the nearest vertex of 2 and vice versa, The 829 Wasserstein-1 distance ( 1 ) or Earth-mover's distance (EMD) is the minimum total area weighted distance 830 of 1-to-1 matching vertices on 1 and 2 . 1 accounts for the area of triangle faces and is minimal if the 831 vertices of 1 is a uniform sampling of 2 or vice versa. Exact computation of 1 is impractical, even for small 832 meshes. We compute 1 using the sliced-Wasserstein, 1 approximation, which uses random spherical given as a percentage. The volume of a mesh was computed as the number of voxels in its binary voxelization 839 computed as described above. We used the minimal possible dilation ball kernel size to ensure a correct 840 volume computation -visual checking of binary voxelization and value at least 3x surface area. We use 1 = 841 topo ( , , ) and 2 = ( , , ) to compute the metrics of Extended Fig. 2,3.

843
Reference surface inference for measurement of protrusion height 844 An optimal reference surface for protrusion segmentation must be a ( , ) parameterized surface i.e. ( ( , , )) for further surface-based processing. We first apply the binary protrusion 896 segmentation above, ( ( , , )) to ( ( , , )), taking the intersection and keeping 897 segmentations with size > 100 voxel 2 Cartesian 3D surface area. We diffuse segmentation labels with 898 labelspreading, clamping ratio 0.99 for 10 iterations, with affinity matrix , = 0.9 as above. We do not 899 rebinarize the label probability at the start of each iteration. Finally, we apply ( ( , , )) to the 900 diffused segmentations to get the final instance segmentation labels, ( ( , , )). (which was not certified by peer review) is the author/funder, who has granted bioRxiv a license to display the preprint in perpetuity. It is The copyright holder for this preprint this version posted April 20, 2023. ; https://doi.org/10.1101/2023.04.12.536640 doi: bioRxiv preprint Equiareal disk parameterization by mesh relaxation. We relax the conformal disk parametrization whilst 917 preserving the boundary topology using the area-preserving flow method 112 . We solve Poisson's equation to 918 compute the smooth vector field for diffusing the area distortion and explicit Euler integration to advect vertex 919 points iteratively with Delaunay triangle flips. The extent of area relaxation achieved is determined by the 920 mesh quality and number of vertices with respect to the extremity of local area distortion. In general, relaxation 921 was less stable compared to our relaxation for spherical surfaces above. For thin and long protrusions, prior 922 downsampling and uniform remeshing of the protrusion submesh was necessary to enable full area distortion 923 relaxation.

925
To convert a unit disk parameterization to an x pixel image, we 'square' the disk using the elliptical grid 926 mapping formula 113 , multiply the resulting vertex coordinates by /2 and interpolate the coordinates and 927 associated vertex quantities onto a x pixel integer grid. This gives similar results to but is significantly 928 faster than solving the Beltrami equation 94 . 929 930

Refinement of undersegmented blebs 931
Given ( , , ) and the vertex ids corresponding to protrusion , ( ) we first impute any small holes in the 932 segmentation; inner vertices not assigned to protrusion but should be in order to ensure the protrusion 128 pixel square image as described above. Positive curvature 'seed' regions are identified by thresholding 941 the mean curvature mapped to 2D, ( protrusion ( , , )) > ℎ ℎ with a global threshold and then applying 942 morphological closing, disk kernel radius 1 pixel. To classify regions as having negative, flat and positive 943 mean curvature, 3-class Otsu thresholding was applied to ( ( , , )) to give two thresholds. All regions 944 with mean curvature greater than the higher threshold ℎ ℎ were positive curvature. Undersegmented 945 blebs correspond to a binary composed of conjoined pseudo-circular regions. We use the gradient 946 watershed 114,115 on the Euclidean distance transform of the high curvature region binary to automatically 947 separate conjoined blebs without seed markers. Mesh matching and interpolation was used to map and 948 segmentation labels between protrusion ( , , ) and protrusion ( , , ) . The refined segmentation were 949 mapped as seed labels from protrusion ( , , ) for every protrusion back to ( , , ). The revised seed labels 950 were then diffused across ( , , ) using the combined geometrical and convexity affinity matrix from above 951 for 10 iterations with = 0.99 . The binary protrusion segmentation from above is applied, and any 952 segmentation with Cartesian 3D surface area < 10 voxels 2 removed to give the final refined protrusion 953 segmentation instances.

972
Volume propagation of surface-based instance protrusion segmentation. The surface-based protrusion 973 segmentation, ( ( , , )) is converted to voxel-based by setting the value of the voxels 974 corresponding to the integer discretized ( , , ) coordinates to the matching protrusion label ID. We expand 975 the labels by 3 voxels using the Python Scikit-Image skimage.segmentation.expand_labels function and 976 mask with the topographic binary cell volume ( ( , , )) to get the initial topographic volume 977 protrusion segmentation, ( ( , , )) with only the surface of protrusions labelled. We apply 978 marker watershed segmentation slice-by-slice to propagate labels laterally into the protrusion volume within 979 a slice and labels from previous slices, from the top, = + to the bottom, = − of the topographic 980 volume. At a slice = 0 , we use the Euclidean distance transform of ( ( = 0 , , )) for watershed 981 with the seed markers given by the labels of the previous slice, = 0 + 1 combined with the current labels 982 of ( ( = 0 , , )) at slice = 0 . In combining labels, the labels of the previous slice = 0 + 1 983 takes precedence and overwrites the label of ( ( = 0 , , )). The result, ( ( , , )) 984 assigns a protrusion ID to all voxels in the entire cell volume ( ( , , )) (Extended Fig. 4f). where on the boundary, we use the mass, boundary and Laplacian, 0 boundary matrix defined for a 2D line 015 and mesh , 0 mesh is the mass and Laplacian matrices defined for a 3D triangle mesh.
[ | ] is used to denote 016 the augmented matrix formed by appending the columns of matrix and . We solve for the vertex position 017 at the next timepoint ( + ) as with cMCF by direct matrix inversion. 018 019 . CC-BY 4.0 International license made available under a (which was not certified by peer review) is the author/funder, who has granted bioRxiv a license to display the preprint in perpetuity. It is The copyright holder for this preprint this version posted April 20, 2023. ; https://doi.org/10.1101/2023.04.12.536640 doi: bioRxiv preprint Surface curvature measurement 020 The mean curvature, of a 3D surface was measured as the divergence of ̂, the unit surface normal 103 , 021 = − 1 2 ∇ ⋅̂. The surface mesh is voxelized to a binary volume, and ̂ is computed as the gradient of the 022 signed distance transform of with the Euclidean distance metric. computed in this manner as opposed 023 to from the mesh directly using discrete differential geometry 103 or quadric plane fitting 117 which is less 024 affected by the number of mesh vertices or the mesh quality. 025 026 Mesh quality measurement 027 The radius ratio = 2 , defined as twice the ratio between inradius and circumradius was used to measure 028 the face quality for a triangle mesh in Extended. Fig. 2,3. It is a mesh quality measure in the sense that the 029 radius ratio obtains its maximum value of 1 for an equilateral triangle; the shape which jointly maximizes all 030 internal angles and gives the best conditioning number for the mesh Laplacian matrix 118 . 031 032

Surface rendering 033
Triangle meshes were exported from Python using the Python Trimesh library into .obj mesh files and 034 visualized in MeshLab 119 . Volumetric images were rendered in Fiji ImageJ through the volume viewer plugin, 035 and intensities were contrast enhanced for inclusion in the figures using Microsoft PowerPoint. The local 036 surface maximum intensity projection image of Fig. 4c was produced by extending z-axis (depth) rays at 037 every xy pixel, and taking the maximum intensity of voxels within ±9 voxels (±1μm) of the cell surface. 038 039 Datasets 040

Cell morphology validation dataset 041
To validate u-Unwrap3D (Fig. 2,3 and Extended Fig. 2-4), we used 66 cell surfaces segmented and surface 042 protrusions classified using u-Shape3D and acquired from high resolution light sheet microscopy 4,120 as 043 previously described 45  We detected ( , ) pixels on blebs by labelling spatially contiguous areas of positive mean curvature based 126 on 3-class Otsu thresholding defining positive, flat or negative curvature. The largest connected component 127 within a bleb bounding box was defined as on-bleb and the remainder area within the bounding box as off-128 bleb. We extracted distortion-corrected average timeseries of bleb area, mean curvature and septin intensity, 129 that is of a scalar quantity, by observing that the mean of over a Cartesian 3D surface area is equivalent 130 to computing a weighted mean over the equivalent ( , ) area, ∬ ( ( , , )) ∬ = ∬ ( ( , )) ∬ . The weight, 131 is the magnitude of the differential area element, = | × | described above. 132 133 Bleb event alignment. 134 Individual blebbing events were detected within a track by applying peak finding after central moving 135 averaging of bleb area timeseries with a window of 3 timepoints. A peak was defined as having a prominence > 136 0.5 and separated from a neighboring peak by at least 3 timepoints. Individual bleb event timeseries were 137 constructed and temporally aligned using the detected timepoint of maximal bleb area as timepoint 0 and 138 taking a window of 14 timepoints on either side (a total 29 timepoints, 35 s). intensities were calculated by extending a trajectory to an absolute depth of 1 μm along the steepest gradient 159 of the distance transform to the mesh surface, and taking the mean intensity along the trajectory. 160 161 u-Unwrap3D analysis. The first timepoint surface mesh was used as the input ( , , ) to u-Unwrap3D to 162 create a common static ( , , ) coordinate space that the surface meshes from all timepoints is mapped to 163 in Step 5 of u-Unwrap3D to generate ( , , ). We use Unwrap-3D with the following parameters for each 164 step: Step 1, cMCF with maximum iterations = 50, = 1 × 10 −5 , automatic stopping threshold, Δ ℎ ℎ = 5 × connectivity and inserting additional triangles to 'stitch' the image boundaries into a spherical topology, (note 178 the latter stitching requires an even number of columns after discounting that the last column is the same as 179 . CC-BY 4.0 International license made available under a (which was not certified by peer review) is the author/funder, who has granted bioRxiv a license to display the preprint in perpetuity. It is The copyright holder for this preprint this version posted April 20, 2023. ; https://doi.org/10.1101/2023.04.12.536640 doi: bioRxiv preprint the first column, hence a 1025 x 512 grid) and = (maximum internal distance transform value ) / steps; 180 for Step 6, marching cubes meshing of the topographic 3D mapped binary cell segmentation at isovalue 0.5 181 following Gaussian smoothing = 1, ACVD remesh with number of clusters = 50% the number of vertices in 182 the marching cubes mesh and retaining the largest connected component mesh. Topographic cMCF with 183 robust mesh Laplacian 102 , mollify factor = 1 × 10 −5 , = 5 × 10 4 was applied to each ( , , ) for 10 184 iterations to compute the corresponding ( , ).

186
Optical flow ROI tracking. We used motion sensing superpixels (MOSES) 83,84 in dense tracking mode which 187 automatically monitors the spatial coverage of ROIs and introduces new ROIs dynamically to ensure uniform 188 spatial tracking at every timepoint. We partitioned the image with an initial user-specified 1000 non-189 overlapping with standard deviation equivalently defined as the square root of the variance, 215 Figure 1. Overview of the surface-guided computing framework u-Unwrap3D. a) Overview of the 6 key 220 steps to map an input genus-X Cartesian 3D surface, ( , , ) such as that obtained from surface meshing 221 the binary segmentation of an input 3D volume image, , and associated scalar measurements, ( ( , , )), 222 via a smooth genus-0 reference surface, ref ( , , ), into any of three additional representations; topographic 223 3D surface, ( , , ), 3D sphere, 2 ( , , ) and 2D plane, ( , ). denotes mesh vertex coordinates, is 224 the mesh mass matrix, is the mesh Laplacian matrix, the step size of area-distortion relaxation, the step 225 size (in pixels) of the propagation distance, the topographic depth (in pixels) and the iteration number. c) 226 Steps to directly remap input genus-0 surfaces without need for a reference surface. In the figure, (⋅), (⋅) 227 denote surface and volume geometries, respectively, in either Cartesian 3D, topographic 3D, or radius-228 standardized 3D spherical coordinates; ( (⋅)) and ( (⋅)) denote surface or volumetric signals as a 229 function of a particular surface or volume geometry. and denote mean and Gaussian curvatures 230 respectively. 231 232 233 . CC-BY 4.0 International license made available under a (which was not certified by peer review) is the author/funder, who has granted bioRxiv a license to display the preprint in perpetuity. It is The copyright holder for this preprint this version posted April 20, 2023. ;

234
. CC-BY 4.0 International license made available under a (which was not certified by peer review) is the author/funder, who has granted bioRxiv a license to display the preprint in perpetuity. It is The copyright holder for this preprint this version posted April 20, 2023. ; https://doi.org/10.1101/2023.04.12.536640 doi: bioRxiv preprint the same cells shown in c). Volume image intensities were visualized using ImageJ volume viewer and 244 contrast-enhanced to better visualize fine protrusions (see Methods). 245 246 . CC-BY 4.0 International license made available under a (which was not certified by peer review) is the author/funder, who has granted bioRxiv a license to display the preprint in perpetuity. It is The copyright holder for this preprint this version posted April 20, 2023. ;

247
. CC-BY 4.0 International license made available under a (which was not certified by peer review) is the author/funder, who has granted bioRxiv a license to display the preprint in perpetuity. It is The copyright holder for this preprint this version posted April 20, 2023. ; https://doi.org/10.1101/2023.04.12.536640 doi: bioRxiv preprint Individual blebs segmented in topographic 3D representation are mapped to Cartesian 3D for visualization 283 and to the 2D plane for tracking (top (which was not certified by peer review) is the author/funder, who has granted bioRxiv a license to display the preprint in perpetuity. It is The copyright holder for this preprint this version posted April 20, 2023. deformation of a closed 3D surface mesh (top) such as by conformalized mean curvature flow (cMCF) 328 (bottom) distorts local geometrical distances and areas as illustrated by zoom-ins. b) Illustration of the two 329 types of metric distortion incurred by mesh deformation; conformal and equiareal. In general, lower conformal 330 error is at the expense of equiareal and vice versa. c) Plot of the conformal error (black dotted line, black left 331 y-axis) and area ratio (red dotted line, red right y-axis) for each iteration of cMCF (Step 1, Fig. 1b) for the 332 example mesh in a) and Fig. 1 with stop iteration indicated by a black vertical dashed line. d) Rendering of 333 the conformal error and area ratio error at each triangle face for quasi-conformal spherical parametrization 334 of the smooth shape to the sphere (Step 2, Fig. 1b). e) Plot of the conformal error (black dotted line, black 335 left y-axis) and area ratio (red dotted line, red right y-axis) for each iteration of the spherical area distortion 336 relaxation with stop iteration indicated by a black vertical dashed line (Step 3, Fig. 1b)  of the total 66 input surface meshes across morphological motif types to u-Unwrap3D (columns 1-6), of which 356 for a total 63 input meshes (>95%) an equiareal spherical parametrization and topographic meshes were 357 successfully computed (columns 7-8, last two columns). In comparison only 6 input meshes (9%) were genus-358 0. An example of a failed mesh is given in c). c) Example of u-Unwrap3D failure when the conformalized 359 mean curvature flow (cMCF) reference shape at the automatic stop iteration has holes too large to be made 360 genus-0 after voxelization and morphological hole closing in our current implementation (see Methods). d) 361 Table summary of the conformal and area distortion error evaluated at the automatic stop iteration of the 362 cMCF (left half, first 4 columns) and the summary mesh statistics of the remeshed cMCF smooth mesh (right 363 half, last 5 columns). Note the area distortion is given as the surface area fraction ratio with the input surface 364 ref ( , , ) as the denominator and equivalent to 1/ used in the area distortion relaxation (see Methods).

365
(see Methods). All tables in b)-g) report numerical values as median±interquartile range. An inf conformal 366 . CC-BY 4.0 International license made available under a (which was not certified by peer review) is the author/funder, who has granted bioRxiv a license to display the preprint in perpetuity. It is The copyright holder for this preprint this version posted April 20, 2023. ; https://doi.org/10.1101/2023.04.12.536640 doi: bioRxiv preprint error indicates local breakdown of flow for a mesh. When this occurs, cMCF cannot continue and the 367 automatic stop iteration is the iteration # just prior to breakdown. 368 . CC-BY 4.0 International license made available under a (which was not certified by peer review) is the author/funder, who has granted bioRxiv a license to display the preprint in perpetuity. It is The copyright holder for this preprint this version posted April 20, 2023.  . 1b). The target and 372 optimal minimum conformal error is 1.0 and is achieved. b) Table summary of the number of iterations, 373 . CC-BY 4.0 International license made available under a (which was not certified by peer review) is the author/funder, who has granted bioRxiv a license to display the preprint in perpetuity. It is The copyright holder for this preprint this version posted April 20, 2023. ; https://doi.org/10.1101/2023.04.12.536640 doi: bioRxiv preprint conformal and area distortion error and mesh quality for four different stopping criteria (labelled i-iv, see 374 Methods) for area-distortion relaxation of 2 ( , , ) → Ω 2 ( , , ). The target and optimal minimum area 375 distortion is 1.0. The implemented spherical area-distortion relaxation scheme in this paper achieves the 376 optimal area distortion, an equiareal parametrization in Ω < a maximum allowed 50 iterations by minimising 377 directly the area-distortion factor, (equiareal) criterion (Methods) (top left, i). The scheme further allows 378 relaxations between conformal and equiareal parametrizations as demonstrated by three additional stopping 379 criteria: with no relaxation, (i.e. iteration 0) the parameterization is conformal and the most isometric 380 parametrization (MIPS) (top right, ii); for iteration numbers ≈ 1 2 Ω , the parameterization minimises jointly 381 the combined conformal and area distortion as measured by + log , the sum of the quasi-conformal 382 error, and the natural logarithm of the area distortion factor, (bottom left, iii); and for iteration numbers 383 ≲ Ω , the parameterization is the area-preserving MIPS (bottom right, iv). c) Table summary of the statistics  384 for computed topographic meshes, ( , , ) using a 1024×512 pixel ( , ) grid for the subset of =63 385 meshes with successful equiareal spherical parametrizations, Ω 2 ( , , ) (Extended Fig. 2b). d) Table  386 summary of the quantitative measurement of geometric error between the Cartesian 3D remapping, 387 topo ( , , ) of the topographic mesh, ( , , ) and the original input mesh ( , , ) for four metrics; ( ) 388 chamfer distance (1 st column), ( 1 ) sliced 1-Wasserstein (2 nd column), Δ , the percentage difference in total 389 surface area (3 rd column) and Δ , the percentage difference in total volume (4 th column). For a perfect 390 reconstruction, all measures should be 0. Units are given as voxels due to heterogeneous pixel resolution 391 amongst input meshes. A large Δ , but small Δ for blebs and filopodia were due to a subset of non- CC-BY 4.0 International license made available under a (which was not certified by peer review) is the author/funder, who has granted bioRxiv a license to display the preprint in perpetuity. It is The copyright holder for this preprint this version posted April 20, 2023. ; https://doi.org/10.1101/2023.04.12.536640 doi: bioRxiv preprint 408 . CC-BY 4.0 International license made available under a (which was not certified by peer review) is the author/funder, who has granted bioRxiv a license to display the preprint in perpetuity. It is The copyright holder for this preprint this version posted April 20, 2023. ; https://doi.org/10.1101/2023.04.12.536640 doi: bioRxiv preprint