Research Article  Open  Published:
The dynamics of growth cone morphology
BMC Biologyvolume 13, Article number: 10 (2015)
Abstract
Background
Normal brain function depends on the development of appropriate patterns of neural connections. A critical role in guiding axons to their targets during neural development is played by neuronal growth cones. These have a complex and rapidly changing morphology; however, a quantitative understanding of this morphology, its dynamics and how these are related to growth cone movement, is lacking.
Results
Here we use eigenshape analysis (principal components analysis in shape space) to uncover the set of five to six basic shape modes that capture the most variance in growth cone form. By analysing how the projections of growth cones onto these principal modes evolve in time, we found that growth cone shape oscillates with a mean period of 30 min. The variability of oscillation periods and strengths between different growth cones was correlated with their forward movement, such that growth cones with strong, fast shape oscillations tended to extend faster. A simple computational model of growth cone shape dynamics based on dynamic microtubule instability was able to reproduce quantitatively both the mean and variance of oscillation periods seen experimentally, suggesting that the principal driver of growth cone shape oscillations may be intrinsic periodicity in cytoskeletal rearrangements.
Conclusions
Intrinsically driven shape oscillations are an important component of growth cone shape dynamics. More generally, eigenshape analysis has the potential to provide new quantitative information about differences in growth cone behaviour in different conditions.
Background
Brain function depends on precisely specified patterns of wiring between neurons, and failures of wiring can compromise normal function [13]. This wiring develops during early life as axons grow and navigate to find their appropriate targets, often over long distances. A critical role in guiding axons to their targets during neural development is played by neuronal growth cones [4,5]. Growth cones have a remarkably complex and dynamic morphology, with their shape changing on the timescale of minutes [68]. In vivo some of these shape changes appear related to the position of the growth cone along its trajectory, with more complex morphology at choice points [914] suggesting that shape changes play an important role in guidance. However, previous morphological analyses of growth cones have been largely driven by human judgement regarding important shape dimensions, rather than these dimensions being determined directly from the data.
The most prominent features of growth cone structure are filopodia and lamellipodia. Filopodia can be quantified in terms of their number, positions, angles and lengths, while a simple measure of lamellipodial extent is the total area of the growth cone. One way of quantifying the shape of a growth cone at each moment is therefore to provide a list of these quantities, which for a typical growth cone with say five filopodia would consist of 21 numbers (two for the position coordinates and one each for angle and length for each filopodium, plus total area). While such a quantification can be useful, it clearly has significant limitations. First, it relies on timelapse imaging of a resolution sufficient to resolve all individual filopodia, which can be difficult to achieve for dynamic growth cones for long periods of time, especially in vivo. Second, despite its length this list still ignores many obvious characteristics of growth cone shape, such as the shape of the lamellipodia. Third, the formulation of this list takes no account of the statistical structure of the actual data: it is driven by human intuition rather than an objective dissection of where the most variance of growth cone shape actually lies.
Here we therefore take a different approach to quantifying growth cone shape. The approach is based on principal component analysis (PCA), a wellknown mathematical method for revealing the dimensions of a dataset which have the most variance. The application of PCA in the context of shapes is often called eigenshape analysis [15]. Each shape can be parameterised by the coordinates of a set of points, usually placed around the perimeter of the shape. This vector of coordinates can be represented as a point in a highdimensional space, so that a set of shapes is represented by a cloud of points in that space. PCA rotates the axes of this space so that the axes are now ordered in terms of the variance in the data they explain. This reveals the directions (dimensions) through this cloud of points along which there is maximum variance, i.e., along which the cloud of points is most spread out. The first few principal components or eigenshapes then represent the most important shape dimensions, which can be seen as the fundamental building blocks of the set of shapes in question. On a cellular scale, eigenshape analysis has previously revealed important information about the shapes of keratocytes [16] and Dictyostelium [17]. It has also proved to be an extremely useful data analysis tool in domains as diverse as Caenorhabditis elegans locomotion [18], computer vision [19], palaeontology [20], botany [21] and musical instrument design [22]. Here we use eigenshape analysis to reveal the basic building blocks of growth cone morphology, previously unknown properties of how growth cone shape evolves through time, and new insights into the relationships between growth cone shape, chemotactic responses and forward movement. We then show that a simple computational model of shape changes based on dynamic microtubule instability can quantitatively reproduce the characteristic timescales present in the data.
Results
Growth cone eigenshapes
To generate a database of growth cone shapes we first made timelapse movies of growth cones from neonatal rat superior cervical ganglion neurites (n=163) growing in vitro for 2 to 8 h (mean 2.6 h) at 15 s to 1 min intervals (see Methods, Table 1 and Figure 1a). From these we determined characteristic growth cone shapes using eigenshape analysis, i.e., PCA in the space of shapes for the dataset [15] (Figure 1b). The outline of each growth cone in each frame (n=25,461) was automatically extracted, and parameterised by 250 evenly spaced points. The vector of 500 numbers (250 coordinate pairs) representing each outline can be represented as a point in a 500dimensional space (approximately 25,000 points in total for this dataset). PCA was then applied to extract the directions in the shape space that captured the largest proportion of the variance. The Bayesian information criterion (BIC) [23], which trades off variance explained versus model complexity (number of principal components retained), can be used to determine objectively the number of these dimensions that capture most of the variance of the set of shapes. According to this criterion, the optimal model retained only the top five components, which captured 86% of the variance in growth cone shape (Figure 1c). We henceforth refer to these as the significant modes.
These shape modes split into three types, which we term reflective (R), symmetric (S) and mixed (M) modes (Figure 1d). For R modes, the shapes corresponding to equal positive and negative movements along the shape axis are approximately reflections of each other, with the principal R mode (R1) representing bending of the growth cone. S modes instead have approximate symmetry about their midline for all positions along the shape axis, with S1 representing spreading of the growth cone. M modes have neither of these properties, but may represent a linear combination of R and S modes (see Figure 1). The split between these different types continues for the higherorder modes (Figure 1e). The number following R, S or M refers to the logical sequence for each type of mode. This is most clear for the R modes: R1 displays one bend, R2 (see later) displays two bends, R3 displays three bends, and so on. The modes were robust to the number of points used to define the outline, provided this exceeded approximately 200 (Figure 2a–d). Thus, just five characteristic shape modes capture most of the variance in this set of growth cones, and these modes describe distinct features of growth cone morphology.
Taken individually, none of these mode shapes explicitly represent filopodia. No information has been lost overall, however: the axes of the space have simply been rotated, and adding together a sufficient number of mode shapes will perfectly reconstruct the original shape to any desired level of accuracy (Figure 1f). Rather, the analysis has determined that choosing a few axes that represent individual filopodia explicitly does not capture the largest amounts of variance in the data. It is easy to see intuitively why this would be the case. To construct a filopodium at a particular location using the PCA axes requires the addition of many higherorder modes, each of which by itself only explains a very small amount of variance across the whole dataset. Individual filopodia are sparse, in the mathematical sense that most of the time there is no filopodium at a particular position, but when it is present it is represented strongly. We therefore also performed an independent component analysis [24], which is well suited for extracting sparse structure in data. As expected, this produced mode shapes hinting at filopodialike structures (Figure 2e). However, these independent component modes are less useful for capturing general patterns of overall shape [15], and we therefore focused on eigenshape analysis.
To maximise throughput, we used a plastic substrate and relatively low magnification imaging for our experiments. Could the significant eigenshapes derived from these data be missing key features of growth cone shape that would become apparent from analysis of higherquality images? To address this we also performed experiments on a glass substrate (ten movies, 2,325 frames), which allowed, for instance, clearer visualisation of filopodia. In this case there were seven significant mode shapes, but their form was similar to those observed on a plastic substrate (Figure 3a). From this dataset we also generated degraded growth cone outlines to match approximately the quality of outlines available using the plastic substrate. The resulting modes shapes were very similar (Figure 3b; in this case there were five significant modes). Thus the leading shape modes are not very sensitive to the level of detail at which the growth cone outlines are captured, making this form of quantification suitable for a much broader range of imaging regimes than quantifications relying on precise identification of individual filopodia.
Growth cones oscillate
Each growth cone outline can be projected onto each shape mode to give a set of mode scores. These measure the degree to which each shape is represented by each outline, i.e., the position of the growth cone shape along that shape axis. The overall mode score frequency distributions are shown in Figure 4. As expected the R mode projections have roughly symmetric frequency distributions, while M modes have slightly asymmetric distributions, and S modes have highly asymmetric distributions. There are no linear relationships between the distributions of pairs of mode scores since, by definition, PCA dimensions are orthogonal. A statistical test for nonlinear relationships [25] showed some dependencies between modes, though plotting pairs of mode scores against each other did not reveal any obvious patterns (data not shown).
We then examined how mode scores varied through time for each individual growth cone (Figure 5a). These revealed strong hints of periodic behaviour. For instance in the top left panel of Figure 5a, the R1 mode score varies regularly between positive and negative values, indicating a periodic alternation in the growth cone shape between bending left and bending right. To detect periodic patterns quantitatively, we used the autocorrelation and Fourier power spectrum of mode scores. Oscillations in mode scores were common in all significant modes (examples shown in Figure 5b,c,d). Similar oscillations were seen on a glass substrate (data not shown). Oscillation strength and frequency for each mode were quantified using a modified version of the method of [26] (see Methods). The distributions of oscillation strengths (hereafter oscillation scores) and frequencies for all significant modes (n=5×163=815) are shown in Figure 6a,b. The mean oscillation frequency was 0.0338 ± 0.0197 min ^{−1} (mean ± standard deviation), corresponding to a mean period of 30 min. The mean oscillation score was 6.0 ± 2.6. R1 oscillations were on average significantly stronger than oscillations in other modes (P<0.03, Wilcoxon ranksum tests). We used shuffled controls to demonstrate that these oscillations were not simply an artefact of our analysis methods (example in Figure 5e, histograms in Figure 6a,b). The mean shuffled frequency was 0.0754 min ^{−1}, and the mean shuffled score was 3.5. ttests comparing the real and shuffled distributions gave P values of 10^{−76} for frequencies and 10^{−107} for scores. Different modes for the same growth cone sometimes oscillated at the same frequency, and sometimes at different frequencies (see examples in Figure 5), with the former case revealing a range of phase relationships (Figure 6c). Overall, we conclude that growth cone shape oscillates on an average timescale of 30 min.
Relationship between oscillations and movement
Given the variability of oscillation strengths and frequencies about their means, we were interested in whether there were any predictive relationships between mode scores, or oscillation strengths and frequencies, and growth cone movement between frames, characterised by step lengths (distance moved) and bearing changes (angular difference in direction) between consecutive frames. Step length from time t to t+1 was weakly positively correlated with S1 mode score at time t (r=0.11, P=10^{−47}, Spearman correlation), meaning that wider growth cones had a small tendency to take larger steps. Other mode scores were almost uncorrelated with the movement of the growth cone (Table 2), and in particular the correlation of the R1 mode score with subsequent bearing change, although significant, was very small. Thus the periodic alternation in shape of bending left versus bending right did not translate into zigzags in overall axon trajectory. Rather, a periodic sweeping left and right of the growth cone was superimposed on a steady forward motion. R1 oscillations could be seen as a way by which growth cones might systematically probe or ‘sniff’ the environment for molecular or topographical cues, sweeping out more area than they would without these oscillations.
However, oscillation strength and frequency for each movie were correlated with the movement of the growth cone averaged over all frames (Table 3). Average growth rate (step length averaged over each movie) was correlated with S1 oscillation strength (see Methods, r=0.20, P=0.016, n=150, Figure 6d) and R1 oscillation frequency (r=0.37, P=10^{−6}, Figure 6e). Performing a linear regression of average step length against the ten variables of oscillation score and frequency for the five significant modes produced a good prediction of average growth rates (R=0.46, P=0.0002, Figure 6f, Table 4). (Note that oscillation strength and frequency are properties of an entire movie, and thus cannot be correlated with growth cone behaviour at one moment in time, such as instantaneous growth rate.) Thus, growth cones that are oscillating strongly and rapidly tend to make forward progress faster.
Oscillations during chemotaxis
An important mechanism by which axons are guided in vivo is chemotaxis. We therefore applied eigenshape analysis to determine the characteristic behaviour of a new set of growth cones as they underwent chemotactic movement in the growth cone turning (or pipette) assay for 1 h [27] (Figure 7a, Table 5). This assay produces a gradient steepness of approximately 10% across the growth cone [28]. Previously it was found that axons growing in gradients of steepness less than 1% show a bellshaped chemotactic sensitivity curve predicted well by a Bayesian model of the chemotactic response [29]. We confirmed that a similar curve holds for steeper gradients (Figure 7b), as theoretically predicted [30]. By reducing PKA activity in the growth cone via addition of KT5720 [3133], we also found that the chemotactic response for repulsion was approximately the mirror image of that for attraction (Figure 7b).
In this dataset there were six significant eigenshape components, explaining 87% of the variance (Figure 7c). These were very similar to those we observed for the original in vitro (no gradient) dataset, but more clearly illustrated the split into R and S modes. These modes oscillated similarly (Figure 7d), though the short assay duration prevented observation of longer oscillations. There were no significant differences between oscillation strength and period between different gradient conditions, and no relationship was observed between final turning angle and oscillations. However, mean mode projection scores (rather than oscillation scores) varied systematically across gradient conditions, demonstrating a direct effect of chemotactic cues on growth cone shape (Figure 8, Table 6). For instance, the mean R1 mode scores roughly followed the same shape as the chemotactic sensitivity curves.
In a further set of turning assay experiments (n=12, data not shown), we imaged growth cones for 1 h without a gradient, then applied a gradient for 1 h, to determine whether the presence of the gradient would change the properties of the oscillations for each growth cone. However, we found no significant differences.
Eigenshapes and oscillations in vivo
To investigate whether similar behaviour is observed in vivo, we performed timelapse imaging, image segmentation and eigenshape analysis as before on mGFPlabelled growth cones (n=27) of zebrafish retinal ganglion cell axons as they navigated across the optic tectum for periods of 1 to 24 h (mean 14.5 h) (Figure 9a) [34]. This revealed six significant shape modes that were very similar to those seen for the in vitro no gradient and pipette datasets (Figure 9b). Oscillations in mode shapes were also qualitatively similar (Figure 9c,d,e). The mean oscillation frequency was 0.0148 ± 0.0134 min ^{−1}. The difference with the in vitro frequency is likely at least partly due to the longer duration of the movies, which allows lower frequency oscillations to be included in the mean (see the simulations in the following section). The mean oscillation score was 6.6 ± 3.2, similar to the in vitro data. Correlations between oscillations and average step lengths were small (data not shown), but this dataset is much smaller than the in vitro dataset. Thus, the eigenshapes and oscillations displayed by growth cones as they navigate in a molecularly and structurally complex in vivo environment are similar to those displayed in a much simpler in vitro environment, arguing that these behaviours represent fundamental characteristics of growth cones rather than an in vitro artefact or the response to a specific environment.
A computational model of oscillations
What events inside the growth cone could be driving the strong periodicity in shape dynamics we have observed? Can we explain both the mean period of 30 min, and the large variability about this average? A critical component of the growth cone cytoskeleton is microtubules [35,36]. These extend from the axon shaft into the body of the growth cone, and sometimes into individual filopodia [37]. Microtubule growth is characterised by dynamic instability, whereby phases of growth are followed by catastrophic collapse, and then a return to the growth phase [38]. Walker et al. [39] constructed a computational model of this phenomenon, with parameter values constrained directly from experimental measurements. Janulevecius et al. [40] then adapted this model to show that the small volumes of cells, and thus the limited supply of free tubulin, could significantly impact on microtubule dynamics. As a minimal model of growth cone shape changes, here we consider two microtubules within a growth cone competing for the same limited supply of tubulin monomers (Figure 10a).
We assume that the extension of one of these microtubules corresponds to a shape deviation in one direction (e.g., bending of the growth cone to the right), while the extension of the other microtubule corresponds to a shape deviation in the opposite direction (e.g., bending to the left). Referring to the lengths of the two microtubules at time t as l _{1}(t) and l _{2}(t), we take the normalised difference in lengths L(t)=(l _{1}(t)−l _{2}(t))/(l _{1}(t)+l _{2}(t)) as an analogue for the projection onto an eigenshape mode. The model and parameters for each microtubule were taken directly from [40] (see Methods), with the growth cone volume assumed to be 30 µm ^{3} (i.e., roughly 6 µm in diameter and 1 µm high). The only interaction between the microtubules was via their competition for the same pool of tubulin monomers.
Growth cones (n=400) were simulated for the equivalent of approximately 2.6 h (the mean length of movies in the in vitro dataset, corresponding to 3×10^{6} events in the Monte Carlo simulation). L(t) was then subsampled to one value every minute of simulated time, the time resolution of most of our in vitro experimental data. As expected, l _{1}(t) and l _{2}(t) tended to compete with each other (Figure 10b). We found that L(t) generally showed periodic oscillations over time, qualitatively resembling the experimental data (Figure 10c,d,e). To quantify this behaviour, oscillation scores and frequencies were calculated for L(t) exactly as for the experimental data, producing distributions over the set of simulated growth cones. Oscillation scores were larger than for the experimental data, with a mean and standard deviation 17.4 ± 6.0, compared to 6.0 ± 2.6 for the experimental data (Figure 10f). This is perhaps not surprising, since there are many additional sources of noise in the data, which are not present in the model. However, surprisingly given the simplicity of the model, there was an extremely close match between the oscillation frequencies of the model and the in vitro data: 0.0333 ± 0.0210 min ^{−1} for the model, compared to 0.0338 ± 0.0197 for the data (P=0.9, ttest, Figure 10g). Thus, the model implicates dynamic microtubule instability as the driving force behind oscillations in growth cone shape.
The model also allows an analysis of how the duration of simulated time affects mean frequencies. Intuitively, longer durations would lead to lower average frequencies, since more distant parts of the underlying frequency distribution are now included in the average. Simulating 100 growth cones for 1.7×10^{7} events (approximately 15 h of real time, similar to the mean for the in vivo dataset) gave a mean frequency of 0.0176 ± 0.0124 min ^{−1}, confirming this intuition. This was not significantly different to the mean oscillation frequency observed in the in vivo data (P=0.12), suggesting that the lower mean frequency observed in our in vivo compared to in vitro data is mainly due to the longer average length of these movies, rather than any fundamental difference in oscillations between in vitro and in vivo settings. This is consistent with the idea that oscillations are intrinsically driven.
Discussion
Neuronal growth cones have one of the most dynamic morphologies of any (sub)cellular system, and this is challenging to quantify. Here we have applied for the first time to growth cones a shape analysis technique that has proved useful in many other contexts, and shown that just a few basic shape primitives capture the vast majority of the variance in growth cone shape. The form of the leading modes themselves is quite intuitive: R1 represents bending while S1 represents thinness versus fatness, whether these are a result of lamellipodial or filopodial outgrowth. By reducing shape to just the list of numbers representing the projections onto the significant modes, it becomes possible to identify patterns in the data that are not otherwise apparent. In particular, by examining how these projections evolve over time, we found that shape oscillations are a key organising principle of growth cone shape dynamics. The forward movement of growth cones has previously been reported to follow cycles of pausing and growth [41,42], and zigzag behaviour in gradients [43]. However, our morphometric [44] analysis reveals for the first time the shapes representing the most important degrees of variance in growth cone morphology (i.e., eigenshapes), how the projections of these shapes vary through time (i.e., oscillations), and how they are correlated with growth cone movement (i.e., growth cones with strong fast shape oscillations tended to extend faster). We observed the same general patterns of eigenshapes and oscillations in vitro and in vivo, suggesting that oscillations are a surprisingly robust and fundamental aspect of growth cone behaviour.
We also presented a computational model of growth cone shape changes based on dynamic microtubule instability. This model is clearly a highly abstracted version of reality, and is intended as simply a minimal model of variation along a shape axis (e.g., bending left versus bending right). It is therefore remarkable that such a simple mechanism produces an excellent match to not only the mean frequency of oscillations in the data, but also the variance in this distribution. This provides strong support for the hypothesis that the basic driver of the growth cone shape oscillations we have observed is dynamic microtubule instability in the context of competition for a limited supply of tubulin monomers. Clearly in reality many other components will be involved, most notably the actin cytoskeleton, and how these work together to determine oscillations remains to be determined.
We have examined growth cones in a relatively featureless in vitro environment, a chemotactic gradient, and traversing the optic tectum in vivo. Eigenshapes and oscillations are remarkably consistent between these different cases. This argues that these properties are intrinsic and do not represent responses to the local environment, such as a way of navigating around local obstacles, consistent with an explanation in terms of intrinsic periodicity of microtubule growth. The complete trajectories of growth cones in vivo often involve navigating through many stages and choice points, for instance retinal axons growing from the eye to the brain [45,46] or callosal axons finding their targets in the opposite hemisphere [47]. Eigenshape analysis provides a method for reexamining more quantitatively exactly how shape changes over such trajectories, and how it is correlated with navigational function. Whether shape oscillations change their properties over complex trajectories and how such changes might be related to the environmental cues present at each moment remain to be determined.
Our tracing of growth cone outlines was performed on relatively lowmagnification phase contrast (in vitro) or fluorescent (in vivo) images, and thus fine details of some filopodia were inevitably lost. However, we have demonstrated directly that the eigenshapes we found are quite robust to image quality (Figure 3). Filopodia are not ignored in our analysis: for instance the presence of a filopodium in the righthand side of the growth cone will show its effect by increasing the score for mode R1 (cf Figure 1f). Eigenshape analysis complements rather than replaces finescale analysis of filopodial dynamics, since it combines filopodial and lamellipodial structure to emphasise overall shape trends rather than finescale structure. Despite recent advances for fluorescently labelled growth cones [48], fine analysis of filopodia is still limited by the difficulty in obtaining large datasets: fully automated image analysis techniques have difficulty with phasecontrast imaging at this level of detail, and handtracing is prohibitively timeintensive. Unfortunately manual observations are not scalable [49], meaning that there is currently a growing mismatch between our ability to manipulate growth cones and to assess the effect of these manipulations on growth cone behaviour using manual techniques.
Conclusions
Eigenshape analysis of growth cones has the potential to provide a novel quantitative understanding of the differences between growth cones in different conditions. These include the effects of the environment as axons navigate towards their targets during development, differences between initial development and regeneration, and differences between mutants and wild type. Overall this work reveals a new dimension to the understanding of the dynamic morphology of growth cones, and potentially opens up novel directions for research into understanding the biological basis of developmental brain disorders.
Methods
All experiments were approved by The University of Queensland Animal Ethics Committee and were performed according to the National Health and Medical Research Council’s animal code of practice.
In vitro imaging
Neurons from superior cervical ganglia (SCG) from P0P4 Wistar rat pups were prepared as previously described [50]. SCGs were incubated in 0.25% trypsin (Gibco) at 37°C for 20 min and then triturated for 10 min. Cells were plated in OptiMEM solution containing 10 µg/mL mouse laminin (Invitrogen) and 0.3 nM nerve growth factor (2.5S mouse NGF, Biosensis) and incubated overnight at 37°C in 35 mm plastic or glassbottomed Petri dishes. Phase contrast images of growth cones were acquired at 20 × for 2 to 8 h at 1 min intervals with a Zeiss Axio Observer.
Imaging in steep gradients of nerve growth factor
SCG neurons were prepared as above. After overnight growth, steep gradients of NGF were produced as reported previously [27,50]. Briefly, axons were positioned with their growth cones 100 µm away from a glass micropipette containing NGF, and with their direction of growth at 45° to the pipette tip. NGF was expelled at 2 Hz using 3 psi to create gradients with a 10% to 15% change in concentration across 10 µm. Phase contrast images were acquired at 20 × for 1 h at 1 min intervals with a Zeiss Axio Observer.
Zebrafish imaging
Adult Tupfel longfin zebrafish were cared for by the Australian Zebrafish Phenomics Facility. Embryos 2.5 to 5 days post fertilisation carrying the BGUG transgene [51] for sparse mGFP labelling were used in the experiments. Nphenylthiourea (0.003%) kept the embryos transparent. Embryos were mounted in low melting point 1.5% agarose (SeaPlaque) in E3 medium on a glass coverslip and maintained at 28.5°C. Confocal image stacks were taken through the tectum (typically 40 to 60 µm, 1 µm intervals) and digitally flattened for analysis. GFP images at 10 and 1 min intervals were taken using a Zeiss LSM 510, through either a 20 × or 40 × objective. From the final timelapse movies, spans 1 to 20 h long with visible growth cones were chosen for analysis.
Semiautomated outline capture
All image processing was carried out using customised software developed using Matlab (Mathworks). To extract growth cone outlines, a quadrilateral region of interest was manually placed over the growth cone for selected frames of each movie and then automatically interpolated for intervening frames. Frames where growth cones were clearly bifurcating were excluded. Images within the region of interest were filtered using the Matlab imagetexturing transforms stdfilt and entropyfilt. The growth cone outline was then manually optimised in approximately five to ten frames, and these candidate outlines were used to train linear support vector machines. This involved calculating local metrics (pixel intensity, size and distance from main body) pertaining to growth cone features that the user accepted or rejected (e.g., missing filopodia and foreign cell matter removal). The remaining movie frames were then automatically resegmented using the support vector machines as a basis for attaching disconnected growth cone segments and features. All remaining imperfections were added or removed from the binary images via manual postprocessing. These imperfections included particles floating in the media that had temporarily attached themselves to the growth cone, or defects in the underlying substrate (e.g., dark lines) that had become inappropriately incorporated into the growth cone outlines using the automated analysis. Additions to the growth cone occurred when the segmentation process failed to attach whole parts of the growth cone (e.g., when the axon became too thin to be segmented appropriately, or filopodia had clearly been missed). These judgements were made by the user performing the segmentation. Normally the total user interaction time required per movie was less than 1 h.
Eigenshape analysis
All axon outlines were rotated to a common axis, and truncated such that their length (to the most distal tip of the growth cone) was 30 µm (20 µm for zebrafish data). For the pipette assay data, all growth cone outlines were mirrored about the yaxis such that the pipette tip was always located to the right of the growth cone. Coordinates of the growth cone representation used for the PCA (250 evenly spaced points) were determined by piecewise cubicHermite polynomial parametric representations (calculated using pchip in Matlab). PCA was performed using princomp in Matlab. The BIC was used to determine the number of PCA dimensions to retain. The last significant mode was taken as the first mode with negative BIC, given by
where k is the number of modes to keep, N is the number of frames, d is the original dimensionality and V is the variance of each mode.
Growth cone centre trajectories
The location of the growth cone centre (GCC) was estimated using a custom automated heuristic. The regions of the growth cone that had the darkest pixel intensity and that were more distant from the rotation base were picked as candidates for the GCC. If there were multiple candidate regions, the most likely was determined as the maximum of the transfer function $Y = (y^{h_{1}} \times I^{h_{2}})\phantom {\dot {i}\!}$ where y is the distance from the base of the rotated growth cone, I is the mean pixel intensity of the candidate group, and h _{1} and h _{2} were 1 and 2, respectively. The final GCC was calculated as the pixel weighted centroid of the candidate region. Clear anomalies in the automated method were corrected by manual manipulation.
Trajectory analysis
GCC trajectories were smoothed using a moving average filter (step size of 3), discarding the first and last four frames to remove edge effects. Step length was measured as the Euclidean distance between successive time points, with average path length being the sum of the step lengths divided by the total number of time points. Average displacement was calculated as Euclidean distance between the first and last frames divided by the length of the movie. The growth cone bearing was calculated as the angle of each step, with the bearing change being calculated as the clockwise angle between two given time points. Initial angles for the pipette assay were calculated using the pipette location as a reference. For other assays the initial angle was calculated as the direction of movement over the first five frames. Final turning angles were measured as the angle formed between the initial direction and the direction of the GCC in the last frame from the GCC in the first frame.
Statistical analysis of relationships
Unless otherwise stated all pairwise comparisons of data were performed using Wilcoxon ranksum tests, while Kruskal–Wallis tests were used for comparisons of multiple groups. Relationships between individual mode scores versus trajectory statistics were calculated using Spearman correlations. Multiple linear regression analysis was used to correlate linear combinations of the mode scores with the trajectory statistics. When performing correlations and regression using oscillation statistics, only the first 2 h of movies were considered to remove bias based on movie length. The relationship between mode score oscillations to turning angle and average path length was calculated using both multiple linear regression and Spearman correlations.
Time series analysis
A fast Fourier transform was calculated using fft in Matlab by centring the signal (mode score or growth cone area) about the mean and zeropadding the signal length to the next power of 2. The autocorrelation function was calculated using autocorr in Matlab. The extent of the oscillatory behaviour of the autocorrelation function was quantified using an oscillation strength metric adapted from [26]. In brief, the symmetric autocorrelation function was calculated and smoothed using a Gaussian kernel with a variance depending on imaging frequency. The central peak was truncated to the level of the next highest peak, and the fast Fourier transform calculated for the resultant signal. The final oscillation strength was defined as the power of the fast Fourier transform peak for the autocorrelation function divided by the mean power of the remaining frequencies.
Computational model
The model was exactly as described in [40], except that there are two microtubules rather than one. In brief a microtubule of length l(t) grows according to
and shrinks according to
where t is time, [ T] is the concentration of free tubulin, and k _{ a }, k _{ d } and k _{ s } are rate constants. The frequencies of catastrophe f _{ c } (switch from the growth phase to the shrinkage phase) and rescue f _{ r } (switch from the shrinkage phase to the growth phase) are given by
where a _{ c },a _{ r },b _{ c } and b _{ r } are constants. The values of all parameters were determined from experimental data [40], and are given in Table One of [40]. There was no spatial component to the model, and thus no time delay in the change in concentration available to one microtubule after binding or release of a tubulin monomer by the other microtubule. Simulations were performed using an eventbased Monte Carlo approach as described in [40], coded in Matlab. Each simulation of 3×10^{6} events, corresponding to 2.6 h of real time for one growth cone, took about 3 min to run on a 2.8GHz Intel Core i7 iMac.
Abbreviations
 BIC:

Bayesian information criterion
 GCC:

growth cone centre
 GFP:

green fluorescent protein
 M:

mixed
 NGF:

nerve growth factor
 PCA:

principal component analysis
 R:

reflective
 S:

symmetric
 SCG:

superior cervical ganglion
References
 1
Yaron A, Zheng B. Navigating their way to the clinic: emerging roles for axon guidance molecules in neurological disorders and injury. Dev Neurobiol. 2007; 67:1216–31.
 2
Lin L, Lesnick TG, Maraganore DM, Isacson O. Axon guidance and synaptic maintenance: preclinical markers for neurodegenerative disease and therapeutics. Trends Neurosci. 2009; 32:142–9.
 3
Stoeckli ET. What does the developing brain tell us about neural diseases?Eur J Neurosci. 2012; 35:1811–7.
 4
GordonWeeks PR. Neuronal growth cones: Cambridge University Press; 2005. http://www.cambridge.org/au/academic/subjects/lifesciences/cellbiologyanddevelopmentalbiology/neuronalgrowthcones.
 5
Lowery LA, Van Vactor D. The trip of the tip: understanding the growth cone machinery. Nat Rev Mol Cell Biol. 2009; 10:332–43.
 6
Goldberg DJ, Burmeister DW. Stages in axon formation: observations of growth of Aplysia axons in culture using videoenhanced contrastdifferential interference contrast microscopy. J Cell Biol. 1986; 103:1921–1.
 7
Stirling RV, Dunlop SA. The dance of the growth cones–where to next?Trends Neurosci. 1995; 18:111–5.
 8
Betz T, Koch D, Stuhrmann B, Ehrlicher A, Käs J. Statistical analysis of neuronal growth: edge dynamics and the effect of a focused laser on growth cone motility. New J Phys. 2007; 9:426.
 9
Bovolenta P, Mason C. Growth cone morphology varies with position in the developing mouse visual pathway from retina to first targets. J Neurosci. 1987; 7:1447–60.
 10
Kaethner RJ, Stuermer CA. Dynamics of terminal arbor formation and target approach of retinotectal axons in living zebrafish embryos: a timelapse study of single axons. J Neurosci. 1992; 12:3257–71.
 11
Myers PZ, Bastiani MJ. Growth cone dynamics during the migration of an identified commissural growth cone. J Neurosci. 1993; 13:127–43.
 12
Buettner HM, Pittman RN, Ivans JK. A model of neurite extension across regions of nonpermissive substrate: simulations based on experimental measurement of growth cone motility and filopodial dynamics. Dev Biol. 1994; 163:407–22.
 13
Halloran MC, Kalil K. Dynamic behaviors of growth cones extending in the corpus callosum of living brain slices observed with video microscopy. J Neurosci. 1994; 14:2161–77.
 14
Hutson LD, Chien CB. Pathfinding and error correction by retinal axons: the role of astray/robo2. Neuron. 2002; 33:205–17.
 15
Pincus Z, Theriot JA. Comparison of quantitative methods for cellshape analysis. J Microsc. 2007; 227:140–56.
 16
Keren K, Pincus Z, Allen GM, Barnhart EL, Marriott G, Mogilner A, et al. Mechanism of shape determination in motile cells. Nature. 2008; 453:475–80.
 17
Tweedy L, Meier B, Stephen J, Heinrich D, Endres RG. Distinct cell shapes determine accurate chemotaxis. Sci Rep. 2013; 3:2606.
 18
Stephens GJ, JohnsonKerner B, Bialek W, Ryu WS. Dimensionality and dynamics in the behavior of C. elegans. PLoS Comput Biol. 2008; 4:e1000028.
 19
Turk M, Pentland A. Face recognition using eigenfaces. In: Proceedings of the IEEE Conference on Computer Vision and Pattern Recognition; 1991. IEEE: 1991. p. 586–91. http://ieeexplore.ieee.org/xpls/abs_all.jsp?arnumber=139758&tag=1.
 20
Lohmann GP. Eigenshape analysis of microfossils: a general morphometric procedure for describing changes in shape. Math Geology. 1983; 15:659–72.
 21
Ray TS. Landmark eigenshape analysis: homologous contours: leaf shape in Syngonium (Araceae). Am J Bot. 1992; 79:69–76.
 22
Chitwood DH. Imitation, genetic lineages, and time influenced the morphological evolution of the violin. PLoS One. 2014; 9:e109229.
 23
Schwartz G. Estimating the dimension of a model. Ann Stat. 1978; 6:461–4.
 24
Hyvärinen A, Oja E. Independent component analysis: algorithms and applications. Neural Netw. 2000; 13:411–30.
 25
Szekely GJ, Rizzo M. Brownian distance covariance. Ann Appl Stat. 2009; 3:1236–65.
 26
Muresan RC, Jurjut OF, Moca VV, Singer W, Nikolic D. The oscillation score: an efficient method for estimating oscillation strength in neuronal activity. J Neurophys. 2007; 99:1333–53.
 27
Lohof AM, Quillan M, Dan Y, Poo Mm. Asymmetric modulation of cytosolic cAMP activity induces growth cone turning. J Neurosci. 1192; 12:1253–61.
 28
Pujic Z, Giacomantonio CE, Unni D, Rosoff WJ, Goodhill GJ. Analysis of the growth cone turning assay for studying axon guidance. J Neurosci Meth. 2008; 170:220–8.
 29
Mortimer D, Feldner J, Vaughan T, Vetter I, Pujic Z, Rosoff WJ, et al. A Bayesian model predicts the response of axons to molecular gradients. Proc Natl Acad Sci USA. 2009; 106:10296–301.
 30
Yuan J, Chan S, Mortimer D, Nguyen H, Goodhill GJ. Optimality and saturation in axonal chemotaxis. Neural Comput. 2013; 25:833–53.
 31
Song HJ, Ming GL, Poo Mm. cAMPinduced switching in turning direction of nerve growth cones. Nature. 1997; 388:275–9.
 32
Forbes EM, Thompson AW, Yuan J, Goodhill GJ. Calcium and cAMP levels interact to determine attraction versus repulsion in axon guidance. Neuron. 2012; 74:490–503.
 33
Sutherland DJ, Pujic Z, Goodhill GJ. Calcium signaling in axon guidance. Trends Neurosci. 2014; 37:424–32.
 34
Simpson HD, Kita EM, Scott EK, Goodhill GJ. A quantitative analysis of branching, growth cone turning and directed growth in zebrafish retinotectal axon guidance. J Comp Neurol. 2013; 521:1409–29.
 35
Tanaka E, Sabry J. Making the connection: cytoskeletal rearrangements during growth cone guidance. Cell. 1995; 83:171–6.
 36
Tanaka E, Ho T, Kirschner MW. The role of microtubule dynamics in growth cone motility and axonal growth. J Cell Biol. 1995; 128:139–55.
 37
Geraldo S, GordonWeeks PR. Cytoskeletal dynamics in growthcone steering. J Cell Sci. 2009; 122:3595–604.
 38
Howard J, Hyman A. Dynamics and mechanics of the microtubule plus end. Nature. 2003; 422:753–8.
 39
Walker RA, O’Brien ET, Pryer NK, Soboeiro ME, Voter WA, Erickson HP, et al. Dynamic instability of individual microtubules analyzed by video light microscopy: rate constants and transition frequencies. J Cell Biol. 1988; 107:1437–48.
 40
Janulevicius A, van Pelt J, van Ooyen A. Compartment volume influences microtubule dynamic instability: a model study. Biophys J. 2006; 90:788–98.
 41
Godement P, Wang LC, Mason CA. Retinal axon divergence in the optic chiasm: dynamics of growth cone behavior at the midline. J Neurosci. 1994; 14:7024–39.
 42
Odde DJ, Tanaka EM, Hawkins SS, Buettner HM. Stochastic dynamics of the nerve growth cone and its microtubules during neurite outgrowth. Biotechnol Bioeng. 1996; 50:452–61.
 43
Ming GL, Wong ST, Henley J, Yuan XB, Song HJ, Spitzer NC, et al. Adaptation in the chemotactic guidance of nerve growth cones. Nature. 2002; 417:411–8.
 44
Zelditch ML, Swiderski DL, Sheers HD. Geometric morphometrics for biologists: Academic Press; 2012. https://www.elsevier.com/books/geometricmorphometricsforbiologists/zelditch/9780123869036.
 45
Holt CE, Harris WA. Position, guidance, and mapping in the developing visual system. J Neurobiol. 1993; 24:1400–22.
 46
Kita EM, Scott EK, Goodhill GJ. Topographic wiring of the retinotectal connection in zebrafish. Dev Neurobiol. 2015. doi:10.1002/dneu.22256.
 47
Suarez R, Fenlon LR, Marek R, Avitan L, Sah P, Goodhill GJ, et al. Balanced interhemispheric cortical activity is required for correct targeting of the corpus callosum. Neuron. 2014; 82:1289–98.
 48
Costantino S, Kent CB, Godin AG, Kennedy TE, Wiseman PW, Fournier AE, et al. Semiautomated quantification of filopodial dynamics. J Neurosci Methods. 2008; 171:165–73.
 49
Brown AEX, Yemini EI, Grundy LJ, Jucikas T, Schafer WR. A dictionary of behavioral motifs reveals clusters of genes affecting Caenorhabditis elegans locomotion. Proc Natl Acad Sci USA. 2013; 110:791–6.
 50
Thompson AW, Pujic Z, Richards LJ, Goodhill GJ. Cyclic nucleotidedependent switching of mammalian axon guidance depends on gradient steepness. Mol Cell Neurosci. 2011; 47:45–52.
 51
Scott EK, Mason L, Arrenberg AB, Ziv L, Gosse NJ, Xiao T, et al. Targeting neural circuitry in zebrafish using GAL4 enhancer trapping. Nat Methods. 2007; 4:323–6.
Acknowledgements
We thank Robert Kerr and Elizabeth Forbes for their early contributions to this work, Pranesh Padmanabhan for help with model simulations, and Muming Poo for many helpful discussions. We also thank the three anonymous reviewers, whose insightful feedback significantly improved the paper. We gratefully acknowledge funding from the Australian Research Council (Discovery Project grant DP11010180), the Australian National Health and Medical Research Council (Project grants 631532, 1043044 and 1083707), and a University of Queensland postdoctoral fellowship (RAF).
Author information
Additional information
Competing interests
The authors declare that they have no competing interests.
Authors’ contributions
GJG designed the project, performed the computational modelling for Figure 10, and wrote the paper, with help from all authors. RAF was responsible for the image processing and overall data management. RAF and DJS performed the eigenshape analysis, and DJS and BAB performed the statistical analyses of the data. AWT performed the chemotaxis assays. ZP and BS performed the in vitro (no gradient) assays. Zebrafish assays were designed by EKS and EMK, and performed by EMK. RAF, DJS, ZP, BS and EMK performed growth cone tracing. All authors read and approved the final manuscript.
Rights and permissions
About this article
Received
Accepted
Published
DOI
Keywords
 Axon guidance
 Neurite growth
 Neural development
 Eigenshape analysis
 Shape analysis
 Brain morphometry
 Oscillations
 Microtubules