Analysis of conditional colocalization relationships and hierarchies in three-color microscopy images

Multimolecular interactions are the basis of cellular functions. High-resolution and super-resolution fluorescence microscopy allow us to routinely image three or more molecular entities at the same time within a cell (Huang et al., 2009; Sahl et al., 2017; Valm et al., 2017). Such multicolor imaging may provide a wealth of information about the spatiotemporal regulation of multimolecular interactions in their cellular context. However, without computational and quantitative analysis, much of this information remains unextracted (Bolte and Cordelières, 2006; Lagache et al., 2015).

Colocalization analysis is a cornerstone technique in cell biology for the analysis of microscopy images. It provides information on the localization of molecules within various subcellular compartments, shedding light on their function (Bolte and Cordelières, 2006; Dunn et al., 2011). While colocalization analysis may lack the spatial resolution necessary for the de novo identification of molecular interactions, it does allow the interrogation of known molecular interactions in their cellular context. This capability complements biochemical and biophysical techniques, which are often employed to establish interactions, but taken out of their cellular context (Caetano et al., 2015; Helmuth et al., 2010; Lagache et al., 2018; Lagache et al., 2015).

Colocalization analysis generally falls into one of two categories: (i) pixel/voxel-based, such as the Pearson correlation coefficient (Aaron et al., 2018; Adler et al., 2008; Bolte and Cordelières, 2006); and (ii) object-based, where first objects are detected or segmented in the images and then the colocalization between objects is analyzed (Lagache et al., 2015; Mascalchi and Cordelieres, 2019). Depending on the nature of the objects and the biological question, object-based colocalization analysis can then employ intensity (Manders et al., 1993), area/volume overlap (Mascalchi and Cordelieres, 2019), or inter-object distances (usually between object centroids; Boutte et al., 2006; Helmuth et al., 2010; Lachmanovich et al., 2003; Lagache et al., 2018; Pastorek et al., 2016). For the latter, colocalization analysis can be thought of as a spatial point pattern analysis task (Diggle, 2014; Lagache et al., 2015). With the advent of single-molecule localization microscopy, which yields data in the form of point coordinates, there have been many developments in colocalization analysis through point pattern analysis (Andronov et al., 2016; Malkusch et al., 2012; Pageon et al., 2016; Williamson et al., 2020).

In spite of the many developments and flavors of colocalization analysis, the overwhelming majority is applicable to only two-color microscopy images, i.e., they investigate only pairs of molecular entities at a time. The few previous studies involving three-color colocalization were mainly focused on the question of whether the three imaged entities colocalize (Fletcher et al., 2010; Stauffer et al., 2018). This limits the type of information that they may reveal, leading to the underutilization of multicolor microscopy images.

Here, we developed an approach for analyzing the colocalization relationships between three molecular entities, termed “conditional colocalization analysis.” It is based on the spatial point pattern analysis of detected objects in microscopy images, applicable to images of both punctate and non-punctate objects. Going beyond the question of whether colocalization is present, it addresses the question of whether the colocalization between two molecular entities is influenced, positively or negatively, by their colocalization with a third entity. For molecules known to interact, the extracted colocalization relationships provide evidence for the cooperativity or competition between interactions, and for the enhancement or suppression of interactions at subcellular locations, domains, or compartments. We benchmarked and validated the method using combinations of molecules known to interact or not interact with each other. To then highlight the capabilities of conditional colocalization analysis, we applied it to investigate the colocalization relationships between a receptor and its downstream adaptor in the context of functionally relevant plasma membrane locations.

To derive measures that reflect the extent to which one molecular entity influences the colocalization between two other molecular entities, we took an object-based approach (Lagache et al., 2015; Mascalchi and Cordelieres, 2019; Simoncelli et al., 2020). Objects in the images (the labeled molecular entities) were detected or segmented, and then their colocalization was assessed based on their nearest neighbor distances. We took an object-based approach for two main reasons. First, for a basic measure of the extent of colocalization, the nearest neighbor distance allowed us to calculate the fraction of one object type colocalized with another object type (Mascalchi and Cordelieres, 2019). By interpreting colocalization fractions as probabilities, we then extended the analysis to calculate conditional probabilities and related measures that reflected the relationships between object colocalizations. Second, object-based analysis offered robust means to computationally assess the significance of any observed colocalization, a major challenge for colocalization analysis (Aaron et al., 2018; Costes et al., 2004). Specifically, in an object-based analysis it is relatively straightforward to (computationally) abolish the relationships between objects, e.g., by object randomization, while maintaining individual object integrity (Aaron et al., 2018; Fletcher et al., 2010; Lagache et al., 2015).

Throughout this work, we will employ the following terminology to refer to the three sets of objects, the colocalization relationships of which are being analyzed. The objects will be referred to as target (T), reference (R), and condition (C). The question that the analysis aims to answer is: How much is the colocalization of target objects with reference objects influenced by target and/or reference colocalization with condition objects? To answer this question, we employed a four-step analysis process. The below descriptions and derivations are for images containing diffraction-limited, punctate objects, the positions of which are described by their center coordinates. In the next section, we present an extension for non-punctate objects.

Suppose there are N_T, N_R,, and N_C target, reference, and condition objects, respectively (Fig. 1 A). First, we determined whether each target and reference object colocalized with a condition object (Fig. 1, A and B). This divided the objects into condition-positive groups (TwC and RwC; “w” indicates “with/positive”) and condition-negative groups (TnC and RnC; “n” indicates “without/negative”).

Two objects were considered colocalized when the distance between them was below a particular colocalization radius. The colocalization radius defined the spatial scale at which the colocalization relationships were explored. It also accounted for image properties such as localization precision (including imprecision introduced by the molecule labeling strategy [Bates et al., 2007]) and registration shifts between different channels. As will be shown later, performing conditional colocalization analysis using a range of colocalization radii yields more robust results than employing one colocalization radius.

Next, we calculated the fractions of target objects in the different groups colocalized with reference objects in the different groups, also using the colocalization algorithm of Helmuth et al. (2010). As in Step 1, two objects were considered colocalized when the distance between them was below a particular colocalization radius. From this, we obtained five colocalization fractions (f_coloc): f_coloc(TwR) (analysis of all target and reference; Fig. 1 C, top row), f_coloc((TwC)wR) and f_coloc((TnC)wR) (analysis focused on condition-positive and condition-negative target; Fig. 1 D), and f_coloc(Tw(RwC)) and f_coloc(Tw(RnC)) (analysis focused on condition-positive and condition-negative reference; Fig. 1 C).

For each colocalization fraction, we calculated a corresponding coincidental colocalization fraction in the absence of any true relationship between target and reference $(fcolocnullTR)$⁠. For this, we kept the reference objects in place and replaced the target objects with points on a grid within the cell area (Helmuth et al., 2010). Using the same colocalization radius as for the original data, we calculated $fcolocnullTR$ corresponding to each of the colocalization fractions mentioned above. This procedure calculated the expected (coincidental) colocalization given the number and spatial distribution of reference objects, i.e., it considered the spatial context of the colocalization analysis (Helmuth et al., 2010).

The colocalization fractions obtained in Step 2 yielded five conditional colocalization measures:

Measure 1. The overall colocalization probability, p(TwR):

N_TwR = number of target objects colocalized with reference objects, irrespective of either’s relationship to condition objects.

Measures 2 and 3. The influence of target colocalization with the condition on target colocalization with reference, p(TwR|TwC) and p(TwR|TnC):

N_(TwC)wR = number of target objects colocalized with both reference and condition objects. N_(TnC)wR = number of target objects colocalized with reference objects but not with condition objects. N_TwC and N_TnC are from Step 1.

Eq. 6a follows from the definition of conditional probability. Eq. 6b follows from the definition of p(TwR&TwC) for the numerator (Fig. 1 D, left branch) and Eq. 1 for the denominator. Eq. 6c follows from Eq. 4.

Eq. 7a follows from the definition of conditional probability. Eq. 7b follows from the definition of p(TwR&TnC) for the numerator (Fig. 1 D, right branch) and the equivalent of Eq. 1 for TnC for the denominator. Eq. 7c follows from Eq. 5.

The comparison of p(TwR|TwC), p(TwR|TnC) and p(TwR) assessed the extent to which target colocalization with reference was influenced by target colocalization with the condition.

Measures 4 and 5. The influence of reference colocalization with condition on target colocalization with reference (p^rs(Tw(RwC)) and p^rs(Tw(RnC))):

N_Tw(RwC) = number of target objects colocalized with reference objects that are themselves colocalized with condition objects. N_Tw(RnC) = number of target objects colocalized with reference objects that are not colocalized with condition objects.

As is, p(Tw(RwC)) and p(Tw(RnC)) could not be directly compared with p(TwR) because the number of reference objects employed to calculate these different probabilities was different. Specifically:

N_RwC + N_RnC = N_R (Fig. 1 B, bottom),

N_Tw(RwC) + N_Tw(RnC) = N_TwR (Fig. 1 C, bottom row),

p(Tw(RwC)) + p(Tw(RnC)) = p(TwR) (Eqs. 8 and 9 vs. Eq. 3).

These measures could be >1, and thus they no longer indicated probabilities (unlike the original fractions in Eqs. 8 and 9). Nevertheless, rescaling factored out the number of reference objects used for analysis, allowing us to compare p^rs(Tw(RwC)), p^rs(Tw(RnC)) and p(TwR) to assess the extent to which target colocalization with reference was influenced by reference colocalization with the condition.

In summary, Eqs. 3, 6, 7, 10, and 11 (in addition to Eqs. 1 and 2 in Step 1) provided us with useful measures for conditional colocalization analysis. Using $fcolocnullTR$ instead of f_coloc in these equations yielded the coincidental counterparts of these measures in the absence of any true relationship between target and reference (p_nullTR and $pnullTRrs$⁠). Conceptually, these measures explored different facets of the joint probability distribution describing the colocalization of the target, reference, and condition (see Materials and methods, Eqs. 12, 13, 14, 15, 16, 17, 18, 19, 20, 21, 22 and related text) to address the specific question of how much the condition influences target-reference colocalization (which the joint distribution, as is, does not).

The nullTR counterparts of the conditional colocalization measures assessed their coincidental values in the absence of a true relationship between target and reference. Similarly, it was important to assess their coincidental counterparts in the absence of any true influence of the condition on target colocalization with reference. For this, we randomized the locations of the condition objects, while keeping the same number of condition objects and not altering the target or reference in any way (see Materials and methods and Fig. S1, A–F). We performed the condition object randomization 100 times for robustness (Fig. S1 G). For each randomization instance, we repeated Steps 1–3 and calculated the coincidental colocalization measures corresponding to that instance. Averaging over the 100 randomization instances then yielded the final coincidental colocalization measures (p_randC and $prandCrs$⁠), which were compared with the original data colocalization measures. This allowed us to assess the significance of any revealed influence of the condition on target colocalization with reference.

In the above analysis scheme, we assumed that objects were diffraction-limited puncta and defined the distance between them as the distance between their center coordinates. While this covers many colocalization questions, many subcellular structures of interest, especially as condition objects, could be larger. Therefore, we generalized the above scheme to the case of non-punctate objects through three modifications. First, we replaced Gaussian fitting with image segmentation to delineate non-punctate objects. Second, to calculate the coincidental colocalization fraction (“nullTR” in Step 2 above) for non-punctate target objects, we employed 100 randomizations of the target objects instead of a grid to preserve the nature of target objects (Fig. S1 H). Third, we generalized the definition of distance between objects by defining the colocalization distance from object A to object B as the average (Euclidean) distance from every pixel in object A to the nearest pixel in object B (Fig. S2).

With this definition, the colocalization distance for two punctate objects was the Euclidean distance between their positions, as above. For the case of the distance from a punctate object to a non-punctate object, the colocalization distance was the distance from the position of the punctate object to the closest pixel of the non-punctate object. Note that the colocalization distance defined as such was not symmetric: the colocalization distance from A to B was not necessarily the same as the colocalization distance from B to A. This was because the colocalization distance considered the shapes and sizes of non-punctate objects in addition to their relative locations, all of which were relevant for the assessment of colocalization (Bolte and Cordelières, 2006; Lachmanovich et al., 2003; Fig. S2). While other metrics could be employed to assess colocalization between non-punctate objects, such as the fraction of area overlap (Mascalchi and Cordelieres, 2019), an advantage of the colocalization distance as defined here was that it was equally applicable to punctate and non-punctate objects. With these modifications, the conditional colocalization analysis scheme was applicable to punctate and non-punctate objects and their combinations.

To assess the accuracy of our conditional colocalization analysis scheme, we tested it on simulated data with known ground-truth colocalization measures. To make the simulations realistic, we took experimental data cell masks and simulated various colocalization scenarios within them. To localize punctate objects, we selected sub-pixel coordinates within each cell mask. For non-punctate objects, we either used paxillin patch segmentations from our experimental data or generated non-punctate objects of similar size and shape to the paxillin patch segmentations (Fig. S3, A–C).

First, we validated the core colocalization analysis between two sets of objects covering three scenarios: (i) punctate target and reference, (ii) punctate target and non-punctate reference, and (iii) non-punctate target and reference. We employed multiple parameter combinations (described in detail in Materials and methods), each simulated 28 times (scenarios i and ii) or 15 times (scenario iii) using 28 or 15 cell masks, respectively. In scenarios i and ii (Fig. 2, A and B), using a colocalization radius of three pixels (which was compatible with how the simulated data were generated), the calculated p(TwR) for each parameter combination was very close to its ground-truth value, and, as expected, it was independent of N_T and N_R. Also as expected, the coincidental colocalization probability p_nullTR(TwR) depended only on N_R, reflecting the fraction of the cell mask area covered by the colocalization zones around the reference objects. In scenario (iii), using a colocalization radius of three pixels, the calculated p(TwR) and p_nullTR(TwR) were close to their input/expected values, but they were overall slightly underestimated (Fig. 2 C). We reasoned that this was because the employed colocalization radius was too small relative to the size of the simulated target objects. Indeed, rerunning the colocalization analysis using a radius of six pixels yielded more accurate p(TwR) estimates, although with a concomitant increase in p_nullTR(TwR) (Fig. S3 D). These results highlight how the colocalization distance between non-punctate objects reflects both their distance and size/shape (Fig. S2), and thus the colocalization radius must take these object properties into account.

Next, we validated the conditional colocalization analysis of three sets of objects. We tested three scenarios: (i) all objects were punctate, (ii) the target and reference were punctate while the condition was non-punctate, and (iii) all objects were non-punctate. For each scenario, we employed multiple parameter combinations (described in detail in Materials and methods), each simulated 15 times using 15 cell masks. We then performed conditional colocalization analysis using a colocalization radius of three pixels. We compared the calculated p(TwC), p(RwC), p(TwR|TwC) and p(TwR|TnC) to their input values, and the calculated p(TwR), p^rs(Tw(RwC)) and p^rs(Tw(RnC)) to their expected values (Eqs. 23, 24, and 25).

We found that the conditional colocalization measures were estimated accurately for all scenarios, with a slightly worse performance in scenario iii (Fig. 2, D–G). As expected, all conditional colocalization measures were higher than their corresponding nullTR, reflecting the simulated colocalization relationships. All nullTR values were very similar to each other, demonstrating that the rescaling employed to calculate p^rs(Tw(RwC)) and p^rs(Tw(RnC)) (Eqs. 10 and 11) properly factored out the number of reference objects used for analysis, placing all colocalization measures on the same scale for direct comparison. Also as expected, p(TwR) was the “weighted average” of p(TwR|TwC) and p(TwR|TnC), with the weights given by p(TwC) and p(TnC) (Eq. 23). Of note, the randC equivalent of all conditional colocalization measures had a median value very close to p(TwR). This was because, upon condition object randomization, the randC condition-positive and condition-negative groups of target and reference were, in fact, random subsets of target and reference objects, yielding randC conditional colocalization measures that converged to p(TwR). The true (“Data”) conditional colocalization measures were then higher or lower than their corresponding randC values, depending on the specific simulated scenario. All in all, these results demonstrate that our analysis scheme accurately estimates relevant conditional colocalization measures for both punctate and non-punctate objects.

A critical parameter for conditional colocalization analysis is the colocalization radius. As discussed above, for the analysis of non-punctate target and reference objects, it must account for their size. But even for punctate objects, it must accommodate properties of the images, such as localization precision, registration shifts between channels, and the labeling density of molecules and complexes. As a result, and especially in conventional light microscopy, the colocalization radius is in the tens to hundreds of nanometers, which is large compared to molecular dimensions. Therefore, it is important to determine the influence of the colocalization radius on conditional colocalization analysis and to determine strategies for obtaining robust, sensitive, and specific analysis results.

To investigate this, we applied our analysis to molecules known to interact with each other (positive control) and to molecules not known to interact with each other (negative control). For the positive control, we chose three proteins that are part of focal adhesions (FAs), namely (activated) β₁-integrin, (phosphorylated) focal adhesion kinase (pFAK), and paxillin (Geiger et al., 2009; Vicente-Manzanares and Horwitz, 2011). Paxillin was taken as a marker for FAs, which served as condition objects. As pFAK interacts with cell surface receptors other than integrins and is involved in signaling pathways other than FA formation (Marlowe et al., 2016; Murphy et al., 2019; Rigiracciolo et al., 2019; Sulzmaier et al., 2014), we reasoned that conditional colocalization analysis may reveal differential colocalization between β₁-integrin and pFAK at FAs vs. away from FAs. For the negative control, we chose the transferrin receptor (TfR), the adaptor protein TSAd (T cell specific adaptor), and clathrin heavy chain (CHC). There are no reports of association between TfR and TSAd. At the same time, CHC is a marker for clathrin-coated structures (CCSs), which mediate TfR internalization into the cell (Liu et al., 2010; Mayle et al., 2012). Thus, we would expect a certain level of association between TfR and CCSs, making them relevant condition objects.

We performed three-color immunofluorescence (IF) imaging of the two combinations of molecules in endothelial cells (ECs) using total internal reflection fluorescence microscopy (TIRFM). This allowed us to focus on the bottom cell surface and minimize apparent colocalization from signals existing in different z-planes (Figs. 3 A and S4 A). For images with a punctate nature (all but paxillin), we determined object positions with subpixel localization using Gaussian fitting (Aguet et al., 2013; Jaqaman et al., 2008; Fig. 3, B and C; and Fig. S4, B–E). For the paxillin images, we segmented the paxillin patches using intensity-based thresholding (Jaqaman et al., 2016), focusing on objects larger than the diffraction limit (Fig. 3, D and E). The punctate nature of target and reference in these datasets allowed us to investigate the effect of the colocalization radius on our analysis results in the worst-case scenario, where the colocalization radius is 1–2 orders of magnitude larger than the object dimensions.

We then performed conditional colocalization analysis using a range of colocalization radii, namely 1–6 pixels (81–486 nm). Throughout this study, we will present the analysis results in two ways. First, as boxplot graphs of the colocalization measures and their coincidental counterparts (nullTR and/or randC), showing individual cell measurements (Fig. 3, F and G). Second, as −log₁₀(P value) tables, summarizing the results of the statistical tests assessing the significance of the observed colocalization measures (Fig. 3 H). For p(TwR), the listed P value is that for comparing p(TwR) to its corresponding nullTR. If p(TwR) is significantly greater/less than its corresponding nullTR, it indicates significant colocalization/repulsion. For all other measures in the table, three tests were performed to assess significance, namely comparison to nullTR, randC and p(TwR). A measure indicates significantly increased/decreased colocalization relative to overall colocalization if it is significantly higher/lower in all three comparisons (using the Dunn–Sidak correction to control for the total type-I error). Thus, the table lists the least significant P value (the smallest −log₁₀(P value)) among the three tests for each measure. Below, we summarize our main findings from applying conditional colocalization analysis to the positive and negative control datasets.

In the positive control dataset, activated β₁-integrin and pFAK exhibited significant colocalization with the paxillin patches (representing FAs, the condition objects; Fig. 3 F). Conditional colocalization analysis also identified robust colocalization relationships between β₁-integrin and pFAK in the context of FAs over the employed range of colocalization radii (Fig. 3, G–J; to be discussed next). In stark contrast, in the negative control dataset, no colocalization relationships were identified between TfR and TSAd, even though both exhibited significant colocalization with CHC (representing CCSs, the condition objects; Fig. S4, F and G). This demonstrates a high level of specificity in the analysis, up to a colocalization radius of at least 6 pixels (486 nm), which we attribute largely to the statistics-based assessment of coincidental colocalization (nullTR and randC).

Conditional colocalization analysis revealed that β₁-integrin and pFAK exhibited significant colocalization with each other, mainly for the subset of pFAK colocalized with paxillin. When investigating the fraction of β₁-integrin (target) colocalized with pFAK (reference), p^rs(Tw(RwC)) was significantly higher than its coincidental counterparts and p(TwR) at all colocalization radii, indicating that β₁-integrin colocalized significantly more with the paxillin-positive pFAK subset than with pFAK overall (Fig. 3, G–I). Conversely, when investigating the fraction of pFAK (target) colocalized with β₁-integrin (reference), p(TwR|TwC) was significantly higher than its coincidental counterparts and p(TwR) at all colocalization radii, indicating that the paxillin-positive pFAK subset colocalized significantly more with β₁-integrin compared with overall pFAK (Fig. 3 J). In contrast to the strong colocalization observed between β₁-integrin and the paxillin-positive subset of pFAK, the overall colocalization between β₁-integrin and pFAK was much weaker: pFAK colocalization with β₁-integrin was not significant at any colocalization radius (p(TwR) in Fig. 3 J), while β₁-integrin colocalization with pFAK was borderline significant at colocalization radii ≥3 pixels (p(TwR) in Fig. 3 H).

These results demonstrate the sensitivity of conditional colocalization analysis and its ability to reveal differential colocalization for different molecule subsets, reflecting the heterogeneity in their interaction landscape. They also highlight how identifying the proper subset of pFAK—that colocalized with FAs—is critical for capturing the strong colocalization relationship between pFAK and β₁-integrin, as pFAK not at FAs is most likely interacting with other receptors (Sulzmaier et al., 2014).

While the colocalization between β₁-integrin and the paxillin-positive subset of pFAK was significant at all colocalization radii, the level of significance varied with colocalization radius in a systematic manner: initially, it increased with colocalization radius, reaching a peak/plateau at colocalization radius of 3–4 pixels, and then it decreased as the colocalization radius increased further (reflected in the increase and then decrease of the −log₁₀(P value) of p^rs(Tw(RwC)) in Fig. 3 H and of p(TwR|TwC) in Fig. 3 J, bottom). The lower significance levels for smaller colocalization radii were probably due to higher levels of false negatives in the colocalization classification (as observed in Fig. 2 C). For larger colocalization radii, lower significance levels were probably due to higher levels of false positives and overestimation of coincidental colocalization (nullTR and randC). An intermediate colocalization radius of 3–4 pixels offered the best balance between false positives and false negatives for this dataset, providing the highest level of significance.

As the optimal colocalization radius is generally not known a priori, these results suggest that analyzing the data using a range of colocalization radii is the optimal strategy for robust colocalization analysis. True colocalization relationships are expected to appear as significant over a range of radii (not necessarily the full range employed), with a pattern of increasing and then decreasing significance as the radius changes from small to optimal to large. This criterion for identifying meaningful colocalization relationships also safeguards against the spurious identification of colocalization relationships at isolated colocalization radii (i.e., not a continuous range of radii) due to the increased false-positive rate associated with multiple statistical tests (Farcomeni, 2008).

Complementary to the strong colocalization observed between β₁-integrin and paxillin-positive pFAK, our analysis revealed a substantial lack of colocalization between β₁-integrin and the paxillin-negative subset of pFAK (p^rs(Tw(RnC)) in Fig. 3, G–I and p(TwR|TnC) in Fig. 3 J). The colocalization for this subset was less than its coincidental counterparts and p(TwR). The reduction in colocalization was more apparent for larger colocalization radii (Fig. 3, H and J), probably because, at larger colocalization radii, the pFAK objects making the paxillin-negative subset were further away from FAs. These observations further emphasize the benefit of performing conditional colocalization analysis using a range of colocalization radii, allowing for the analysis of different subsets of target and reference objects. They also suggest that consistent, complementary trends for complementary colocalization measures (e.g., increased p(TwR|TwC) is consistent with decreased p(TwR|TnC), and similarly for p^rs(Tw(RwC)) vs. p^rs(Tw(RnC))) enhance the evidence for the revealed colocalization relationships.

The accuracy of our conditional colocalization analysis is expected to depend on the accuracy of the object detections and/or segmentations, which in turn depend on the images' signal-to-noise ratio (SNR). Therefore, we investigated the performance of conditional colocalization analysis for images with different SNR.

We took the positive control dataset, which had an SNR in the range of 4–10 (see Materials and methods for definition of SNR) and added Gaussian white noise to the images in order to reduce the SNR to 4, 3, and 2 (Materials and methods; Fig. 4 A). After detecting and segmenting the images as above, we performed conditional colocalization analysis (using a colocalization radius of three pixels) and compared the results at lower SNR to those of the original images (Fig. 4, B and C). We found that colocalization between β₁-integrin and the paxillin-positive subset of pFAK was captured at all tested SNRs. However, the significance of p^rs(Tw(RwC)) exhibited a decreasing trend as the SNR was reduced. This implies that, for weaker colocalization relationships, the relevant colocalization measures may lose significance at lower SNRs.

The comparisons also revealed that all measures, together with their coincidental counterparts, decreased with decreasing SNR (Fig. 4 B). This was because, at lower SNR, a smaller number of objects was detected in each channel, as only brighter objects were detectable above the noise level (Fig. 4, D and E). Also, an additional conditional colocalization measure, p^rs(Tw(RnC)), became significant at lower SNR (Fig. 4, B and C). This trend was also most likely because of the retention of only brighter objects at lower SNR; performing conditional colocalization analysis on the original data, but taking only the top 75, 50, and 25% brightest objects, revealed the same trend (Fig. 4 F). The significant p^rs(Tw(RnC)) indicated that β₁-integrin was less likely to colocalize with paxillin-negative pFAK than with the overall, complementing the significant preferential colocalization of β₁-integrin with paxillin-positive pFAK, as observed for the original data at colocalization radii ≥4 pixels (Fig. 3 H).

Overall, these results demonstrate that conditional colocalization analysis has the sensitivity to identify colocalization relationships at SNRs as low as 2. Nevertheless, imaging conditions that maximize the completeness and accuracy of detection/segmentation of the objects of interest will ensure maximal sensitivity and accuracy of the conditional colocalization analysis results.

The spatial compartmentalization of molecular interactions is an emerging theme in signaling. For example, cell surface receptors have been observed to preferentially interact with downstream effectors at particular membrane microdomains and to signal differently from different subcellular compartments (Birch et al., 2021; Delos Santos et al., 2015; Eichmann and Simons, 2012; Jaqaman and Ditlev, 2021; Sorkin and von Zastrow, 2009; Sungkaworn et al., 2017). As a result, there is increasing interest in characterizing the spatial variation of colocalization between molecules within the cell (Gorlewicz et al., 2020; Jaskolski et al., 2005; Taylor et al., 2018). Conditional colocalization analysis is ideally suited to shed light on the spatial regulation and/or compartmentalization of molecular interactions within their cellular environment in an automated and quantitative fashion.

To illustrate this capability, we applied conditional colocalization analysis to the receptor tyrosine kinase VEGFR-2 (Vascular Endothelial Growth Factor Receptor 2) and its downstream adaptor TSAd (also known as VEGF receptor adaptor protein [VRAP]), in the context of FAs and CCSs (Figs. 5 and 6). VEGFR-2, expressed on the surface of ECs, is the main receptor for VEGF (short for VEGF-A), where it initiates signaling pathways that are essential for vascular function (Ferrara et al., 2003; Simons et al., 2016). Its downstream adaptor TSAd mediates a signaling cascade that modulates FA turnover and cell adhesion, critically important for EC motility and angiogenic sprouting (Abedi and Zachary, 1997; Birukova et al., 2009; Gordon et al., 2016; Sun et al., 2012; Wu et al., 2000). Interestingly, IF imaging of VEGFR-2 and paxillin (as a marker for FAs) indicated that a subset of VEGFR-2 tends to localize at/near FAs (Fig. 5 A), raising the possibility that this facilitates VEGFR-2–TSAd interactions to regulate FAs. As for CCSs, VEGFR-2 is internalized from the cell surface via clathrin-mediated endocytosis, in both ligand-independent and ligand-dependent manners, where it can continue downstream signaling (Basagiannis et al., 2016; Lampugnani et al., 2006; Lanahan et al., 2010; Simons et al., 2016). Thus, VEGFR-2 is also expected to colocalize with CCSs, potentially together with TSAd. Conditional colocalization analysis allowed us to investigate the colocalization relationships between VEGFR-2 and TSAd in the context of these two functionally-relevant membrane domains.

We performed three-color IF imaging of VEGFR-2, TSAd, and paxillin or VEGFR-2, TSAd, and CHC using TIRFM, in the absence or presence of 2 nM VEGF (5 min stimulation before fixation, enough for VEGFR-2 activation [da Rocha-Azevedo et al., 2020]) (Figs. 5 A and 6 A). As above, punctate signals (all but paxillin) were localized using Gaussian fitting (Aguet et al., 2013; Fig. 5, B and C; and Fig. 6, B–E), while paxillin patches were segmented using intensity-based thresholding (Jaqaman et al., 2016; Fig. 5, D and E). We then performed conditional colocalization analysis using a range of colocalization radii (1–6 pixels, as above).

First, we measured the extent of VEGFR-2 and TSAd colocalization with FAs, and whether they depended on each other (in a statistical sense) for FA colocalization. Taking VEGFR-2 as the target, FAs as the reference, and TSAd as the condition, only p(TwR) was significant (Fig. 5 F; median value of p(TwR) using a colocalization radius of 3 pixels = 0.27). This was at all colocalization radii, in both the absence and presence of VEGF. This indicated strong colocalization of VEGFR-2 with FAs, regardless of whether VEGFR-2 was colocalized with TSAd, and regardless of stimulation with VEGF. In contrast, taking TSAd as target and VEGFR-2 as the condition (while keeping FAs as reference), we found that p(TwR) was largely not significant, while p(TwR|TwC) was significantly greater than its counterparts for a range of radii, and p(TwR|TnC) was significantly lower than its counterparts for larger radii (as observed for the positive control; Fig. 5, G and H). In addition, p(TwR|TwC) had more significant P values and for a larger range of colocalization radii in the presence of VEGF vs. its absence (Fig. 5 H). These results indicated that only the VEGFR-2-positive subset of TSAd colocalized with FAs, with an increase in colocalization upon stimulation with VEGF. They suggest that VEGFR-2 plays a role in recruiting and/or retaining TSAd at FAs, and that this is enhanced upon VEGFR-2 activation (which promotes VEGFR-2–TSAd interactions).

What about VEGFR-2-TSAd colocalization? Is it influenced by location relative to FAs? For this, we used FAs now as the condition, and then once VEGFR-2 as target and TSAd as reference, and once the converse. In both cases, VEGFR-2-TSAd colocalization was only significant for the subset of TSAd colocalized with FAs (p^rs(Tw(RwC)) in Fig. 5 I and p(TwR|TwC) in Fig. 5 J). The colocalization between VEGFR-2 and the FA-positive subset of TSAd was stronger upon VEGF stimulation than in its absence (more significant P values for a larger range of colocalization radii). These results suggest that not only does VEGFR-2 play a role in TSAd colocalization with FAs, but also FAs play a role in TSAd–VEGFR-2 colocalization. In other words, there appears to be cooperativity in the colocalization of TSAd with VEGFR-2 and FAs.

As for CCSs, conditional colocalization analysis revealed a different interplay between them and VEGFR-2–TSAd colocalization. First, in terms of colocalization with CCSs, VEGFR-2 and TSAd exhibited (statistical) interdependencies almost the opposite of those observed with FAs; TSAd exhibited overall significant colocalization with CCSs, regardless of its colocalization status with VEGFR-2 (only p(TwR) is significant in Fig. 6 F), while VEGFR-2 colocalization with CCSs was enhanced when it was also colocalized with TSAd (p(TwR) and p(TwR|TwC) are significant in Fig. 6, G and H). Second, in terms of their colocalization with each other, VEGFR-2–TSAd colocalization was significant only for the subset of VEGFR-2 colocalized with CCSs (p(TwR|TwC) in Fig. 6 I and p^rs(Tw(RwC)) in Fig. 6 J), which was opposite of the case of FAs. Third, while VEGFR-2-TSAd colocalization increased upon stimulation with VEGF (as indicated by more significant P values for a larger range of colocalization radii; Fig. 6, I and J; consistent with observations in context of FAs), TSAd and VEGFR-2 colocalization with CCSs was reduced upon VEGF stimulation (as indicated by overall less significant P values; Fig. 6, F and H). The reduced colocalization with CCSs upon VEGF stimulation suggests that many activated VEGFR-2–TSAd complexes are already internalized after 5 min of stimulation (Lampugnani et al., 2006; Lanahan et al., 2010) and/or CCSs are not the primary mechanism for VEGFR-2 internalization upon stimulation (Basagiannis et al., 2016).

Altogether, these analyses indicate that VEGFR-2-TSAd colocalization varies by location at the plasma membrane. It is enhanced at FAs and CCS. Yet the two membrane domains exhibit differential interplay with VEGFR-2-TSAd colocalization. In the case of FAs, VEGFR-2 and FAs mutually enhance each other’s colocalization with TSAd. In the case of CCSs, TSAd and CCSs mutually enhance each other’s colocalization with VEGFR-2. To the extent that colocalization is reflective of interactions between VEGFR-2 and TSAd, this raises the question of what leads to this differential interplay of colocalization. One possibility is the presence of different interaction partners or feedback from downstream signaling that modulate these interactions at different membrane domains. Performing further conditional colocalization analysis with other interaction partners and combining it with perturbation of the different molecular players is expected to help determine the molecular mechanisms behind the observed increase in VEGFR-2-TSAd interactions at FAs and CCSs. Importantly, the increased interactions most likely facilitate VEGFR-2 signaling in the right place at the right time to regulate the cellular response to VEGF.

We have developed a new approach, conditional colocalization analysis, to assess the colocalization relationships and hierarchies between three molecular entities from multicolor microscopy images. Going beyond the question of whether two entities colocalize, it addresses the question of whether the colocalization between two entities is increased or decreased by a third entity. It exhibits high specificity and sensitivity, especially when integrating analysis results over a range of colocalization radii and complementary conditional colocalization measures. The source code for performing conditional colocalization analysis is available at https://github.com/kjaqaman/conditionalColoc.

Conditional colocalization analysis of multicolor images offers an approach to shed light on the dependencies, hierarchies, and regulation of known molecular interactions within their native, unperturbed environment. It complements the classic approach of system perturbation to probe these questions in at least two ways. First, because it extracts information from the unperturbed system, it can provide unique insights in cases where molecular perturbations have side effects that complicate the interpretation of results, or where there is functional redundancy between molecules (Wagner, 2000; Zhang, 2012). Second, it enables the initial dissection of a system to gather evidence for the important players in regulating the interactions of interest, thus guiding molecular perturbation experiments. Conditional colocalization analysis can be tailored to particular biological questions through question-specific selection of the objects analyzed, e.g., objects from a particular subcellular region or objects of a certain size or shape.

While our applications of conditional colocalization analysis have been focused on the molecular scale, this analytical approach is applicable to any biological system. Non-punctate objects can take any form or shape (Fig. S5), as long as they are relatively discrete and are segmented properly, e.g., using ridge detection for curvilinear objects (Kittisopikul et al., 2020) or using deep learning approaches (Lucas et al., 2021). Punctate objects can be at a relatively high density in the image, in which case their center coordinates can be determined using advanced computer vision approaches designed to handle high-density point emitters (Nehme et al., 2018; Zhu et al., 2012). Accurate detection and/or segmentation of the objects to be analyzed increases the accuracy and sensitivity of conditional colocalization analysis.

Furthermore, conditional colocalization analysis can be applied to super-resolution microscopy images, which brings the employed range of colocalization radii closer to molecular dimensions (Huang et al., 2009; Sahl et al., 2017). This will most likely deepen the insight that conditional colocalization analysis provides about molecular interactions and their regulation. Of note, for the application of conditional colocalization analysis to super-resolution data of the localization microscopy type, fluorophore blinking must be corrected for a priori, in order to avoid molecule overcounting and artifactual clustering (Hummer et al., 2016; Saguy et al., 2021; Williamson et al., 2020).

Conditional colocalization analysis can be extended in multiple ways. (1) The same analysis can be applied to 2D live-cell images to determine the temporal evolution of colocalization relationships. (2) The analysis can be readily extended to 3D, as long as there is sufficient axial resolution. (3) The analysis can be extended to investigate the colocalization relationships and hierarchies between more than three entities, as long as informative conditional probabilities and related colocalization measures can be defined. As extending the analysis to more than three colors will increase the total number of possible object combinations, the analysis procedure in this case may be iterative to help identify biologically relevant object combinations.

All in all, conditional colocalization analysis is a broadly applicable approach for the analysis of multicolor microscopy images. By quantifying the colocalization relationships between multiple molecular entities, it aids with the dissection of the complex landscape of multimolecular interactions in the cell.

Human telomerase (hTERT)-immortalized microvascular endothelial cells isolated from human foreskin (TIME cells; ATCC) were grown in ATCC’s vascular cell basal medium supplemented with microvascular endothelial cell growth kit-VEGF and 12.5 μg/ml blasticidine (Sigma-Aldrich) for 48 h at 37°C + 5% CO₂ until reaching 80–90% confluence. At confluence, cells were passaged and 5.4 × 10⁴ cells were plated on fibronectin-coated (10 µg/ml; MilliporeSigma), base/acid cleaned, 0.17 mm (#1.5) glass-bottom dishes (14 mm glass diameter; MatTek) for 18 h prior to experiments.

Plated cells were washed once in wash buffer (HBSS + 1 mM HEPES and 0.1% NGS) to remove the culture medium and then incubated with HBSS containing or not 2 nM VEGF-A₁₆₅ (Genscript) for 5 min at 37°C. The incubation solution was then removed and cells were fixed with a 3.2% paraformaldehyde solution made in PBS (Electron Microscopy Sciences) for 15 min at RT. Samples were then washed three times (5 min each) at RT with wash buffer and then permeabilized for 1 min with cold 0.01% Triton X-100 solution made in PBS. After three washes, samples were blocked for 15 min in blocking buffer (1% BSA, 5% NGS in wash buffer) and then incubated for 1 hr with primary antibodies (Table 1) at RT. After three washes, samples were incubated with secondary antibodies (Table 1) for 15 min at RT. After three more washes, they were incubated with AlexaFluor647-conjugated primary antibodies (Table 1) for 15 min at RT. Finally, after three subsequent washes, dishes were incubated with imaging buffer (Oxyfluor 1%, Glucose 0.45%, Trolox 2 nM) to reduce photobleaching.

Cells were imaged at 37°C using an Olympus IX83 TIRF microscope equipped with a Z-Drift Compensator and a UAPO 100X/1.49 NA oil-immersion TIRF objective (Olympus). The microscope was equipped with an iXon 888 1k × 1k EMCCD Camera (Andor; Oxford Instruments). With an additional 1.6× magnification in place, the pixel size in the recorded image was 81 nm × 81 nm. Using the Olympus CellSens software, excitation light of 640, 561, and 491 nm from an Olympus CellTIRF-4Line laser system was directed to the sample at a penetration depth of 90 nm by a TRF8001-OL3 Quad-band dichroic mirror. Fluorescence of different wavelengths was collected, filtered with emission filters of ET520/40m, ET605/52m, and ET705/72m (Chroma), and projected onto different sections of the camera chip by an OptoSplit III 3-channel image splitter (Cairn Research). The different channels were excited and recorded sequentially in the order 640, then 561, and then 491. Images were acquired with MetaMorph (Molecular Devices). Camera EM gain was set to 100 for all acquisitions.

For all datasets, except for activated β₁-integrin/pFAK/paxillin, each channel was imaged once using an exposure time of 99 ms and a power of 3.6, 6.4, and 5.9 mW at the coverslip (measured with a Thorlabs power meter at an incident angle of 0°) for the 491, 561, and 640 nm lasers, respectively. For the activated β₁-integrin/pFAK/paxillin combination, two types of images were collected, before and after photobleaching. Specifically, β₁-integrin and pFAK were imaged in a streaming mode for 300 frames (at 10 Hz) with laser power increased to 5.4 mW for the 491 channel and 9.6 mW for the 561 channel to photobleach the signal. The first frame of each stream, where the β₁-integrin and pFAK signals were largely non-punctate, was used for conditional colocalization analysis of all non-punctate objects (Fig. S5). To also obtain punctate images for β₁-integrin and pFAK, a single image was acquired for each of them after photobleaching, and this was used for colocalization analysis of punctate target and reference objects together with non-punctate condition objects (Figs. 3 and 4). Paxillin (the 640 channel) was imaged only once, as with the other molecule combinations.

For every three-channel image or movie, a brightfield snapshot of the imaged cell region was also acquired, in order to aid with manual delineation of the region of interest (ROI) mask for the ensuing analysis.

All image and data analysis and simulation tasks were performed in Matlab R2020a (The MathWorks). All employed code and software packages are compatible with Matlab R2020a on a Linux 64-bit operating system. Images were loaded into Matlab using Bio-Formats (Linkert et al., 2010).

Punctate objects were detected using the “point-source detection” particle detection algorithm in u-track (https://github.com/DanuserLab/u-track; Aguet et al., 2013; Jaqaman et al., 2008). In brief, the algorithm consists of two steps: (i) a filtering step to determine pixels likely to contain objects and (ii) a Gaussian fitting step to determine the object positions with sub-pixel localization. With the appropriate, wavelength-dependent standard deviation, a two-dimensional Gaussian is a good approximation of the microscope's point spread function (Thomann et al., 2002; Zhang et al., 2007). Largely default parameter values were used, except for the α-value for determining object detection significance by comparing the fitted Gaussian amplitude to the local background noise distribution (Table 2). The α-value was chosen based on the visual assessment of the detection results with the goal of minimizing both false positives (superfluous detections) and false negatives (missed particles). The number of objects detected within the ROI for each dataset is listed in Table 3.

Non-punctate objects were segmented using an intensity threshold to separate foreground from the background (Jaqaman et al., 2016). The threshold was determined for each image individually. For paxillin, it was taken as the 90th percentile of each image’s intensity distribution after noise filtering (using a Gaussian kernel with standard deviation = 1 pixel) and local background subtraction. For β₁-integrin and pFAK (when the non-punctate version of their images was analyzed), the threshold was taken as the 85th percentile, also after the above-described pre-processing steps. The local background was estimated for each image by filtering it with a Gaussian kernel with standard deviation = 10 pixels. The threshold was then applied to the noise-filtered and background-subtracted image to segment the patches. Patches with an area of at least 30 pixels (in the case of paxillin and pFAK) or 10 pixels (in the case of β₁-integrin) were retained for further analysis. The number of segmented patches within the ROI and their areas are listed in Tables 3 and 4. The slight differences in paxillin segmentation properties for the β₁-Integrin/pFAK/paxillin dataset between the two tables are because of the slightly different ROIs employed for analysis when all channels of this dataset were non-punctate (Table 4) vs. when only the paxillin channel was non-punctate (Table 3).

ROI masks, covering the imaged area of each cell of interest, were delineated manually. For TfR/TSAd/CHC and VEGFR-2/TSAd/CHC, delineation was based on the brightfield image. For VEGFR-2/TSAd/paxillin and β₁-integrin/pFAK/paxillin, delineation was based on the brightfield image as well as the paxillin and pFAK images (where acquired).

Our implementation closely followed the original algorithm proposed in Helmuth et al., (2010). A target object was considered to be colocalized with a reference object if the colocalization distance between it and its nearest neighbor reference object (see Results section “Generalization of conditional colocalization analysis to non-punctate objects” for colocalization distance definition) was smaller than a user-defined threshold (“colocalization radius”). The colocalization radius accounted for various considerations, such as registration shifts between channels, segmentation accuracy, average object size, etc. For a given subset of target objects, the fraction of target colocalized with reference (f_coloc) was then calculated as the ratio of the number of target objects colocalized with reference objects to the total number of target objects in that subset.

To calculate the coincidental colocalization fraction of the target with the reference (⁠$fcolocnullTR$⁠) given the cellular context, in particular the number, spatial distribution, size, and shape of reference objects within the ROI mask, in the case of punctate target objects, we followed Helmuth et al. (2010) and replaced the target positions with a grid of positions covering the whole ROI. Specifically, we used each pixel within the ROI as one grid point (i.e., target position), and repeated the above distance calculations, classification, and fraction calculation. In the case of non-punctate target objects, using a grid was no longer applicable because it would not preserve the nature (i.e., general size and shape) of target objects. Therefore, instead of a grid, in this case, we randomized the target object positions within the ROI mask (see next section for details) and repeated the above distance calculations, classification, and fraction calculation. This was performed 100 times to minimize fluctuations (Fig. S1 H).

Objects were randomized within the ROI mask for two purposes. (1) To calculate the conditional colocalization probabilities expected by chance, i.e., in the absence of any true influence of the condition objects on target-reference colocalization, the condition objects (whether punctate or non-punctate) were randomized within the ROI mask. (2) To calculate the coincidental colocalization fraction of target with reference in the case of non-punctate target, the target objects were randomized within the cell mask (as mentioned above). For both cases, randomization was performed 100 times to minimize fluctuations (Fig. S1, G and H).

For randomizing punctate condition objects, to ensure that none of the randomly chosen locations were too close to the ROI mask boundary, the ROI mask was first eroded using a square structuring element with a size equal to the colocalization radius plus three pixels. The new condition object locations were then chosen within the eroded ROI mask, using the same number of condition objects as in the original data.

For randomizing non-punctate condition objects or non-punctate target objects, we conserved not only the number of objects, but also the shape of each object. To optimize the selection of random locations, objects in this case were placed in random locations sequentially, starting with the largest object and going down to the smallest object. The idea behind this order was that it was easier to place small objects between large objects than the other way around. Once an object was placed, the pixels it occupied were removed from the ROI mask, as those pixels were no longer available for any other object (Fig. S1, A–F). In addition, before placing each object, the mask was eroded using a disk-structing element with a radius equal to the average radius of the object being placed (Fig. S1, A–F). The purpose of this was to avoid placing any new object too close to already-placed objects (leading to overlap) or to the ROI mask boundary (leading to part of the object being outside the ROI mask). With this strategy, we were able to randomize the locations of non-punctate objects within the ROI mask while avoiding any overlap between them and while avoiding the scenario of not enough contiguous space to place larger objects (a scenario that happened often when smaller objects were placed before larger objects).

The validation simulations were performed within cell ROI masks taken from our experimental data to keep the simulations as realistic as possible. Several masks were used as repeats for each parameter combination, as described in Results. Cell masks were eroded using a square element with side length = 10 pixels to avoid molecules being positioned too close to the cell edge. Below, we list the parameter combinations employed for validation testing and the simulation strategies for the various scenarios.

(i) Punctate target and reference. N_T = 150, 200, 250, 300, 350, N_R = 150, 200, 250, 300, 350 and p(TwR) = 0, 0.25, 0.5, 0.75, 1. Total: 125 parameter combinations. (ii) Punctate target and non-punctate reference. N_T = 150, 200, 250, 300, 350 and p(TwR) = 0, 0.25, 0.5, 0.75, 1. Reference objects were taken from our experimental data paxillin patch segmentations. Total: 25 parameter combinations. (iii) Non-punctate target and reference. Reference objects were taken from our experimental data paxillin patch segmentations and combined with ellipsoidal target objects that mimicked paxillin patches (Fig. S3, A–C). N_T = 85 (lower than the N_T employed above to accommodate the larger size of non-punctate target objects) and p(TwR) = 0, 0.25, 0.5, 0.75, 1. Total: 5 parameter combinations.

(i) Punctate target, reference, and condition. Given the observation that the estimated p(TwR) was equally accurate for all tested N_T and N_R (Fig. 2, A and B), N_T and N_R were fixed at 250 (middle of the range tested for two sets of objects). Otherwise, N_C = 250, 300 and 350, p(TwC) = 0.25, 0.5, 0.75, p(RwC) = 0.25, 0.5, 0.75, p(TwR|TwC) = 0.25, 0.5, 0.75, and p(TwR|TnC) = 0.25, 0.5, 0.75. Total: 243 parameter combinations. (ii) Punctate target and reference, and non-punctate condition. N_T and N_R were fixed at 250. Condition objects were taken from our experimental data paxillin patch segmentations. p(TwC) = 0.25, 0.5, 0.75, p(RwC) = 0.25, 0.5, 0.75, p(TwR|TwC) = 0.25, 0.5, 0.75 and p(TwR|TnC) = 0.25, 0.5, 0.75. Total: 81 parameter combinations. (iii) Non-punctate target, reference, and condition. Condition objects were taken from our experimental data paxillin patch segmentations and combined with N_T = 85 non-punctate, ellipsoidal target objects and N_R = 90 non-punctate, ellipsoidal reference objects (Fig. S3, A and B). p(TwC) = 0.25, p(RwC) = 0.25, and two combinations of p(TwR|TwC) and p(TwR|TnC), namely 0.75 and 0.25 and the reverse. Total: 2 parameter combinations.

(i) Punctate target and non-punctate reference. To simulate non-punctate reference objects, the paxillin patch segmentations belonging to the employed cell masks were used as reference objects. Then, based on the input p(TwR), target objects were divided into colocalized or not colocalized with reference objects. Target objects colocalized with reference objects were placed 0–2 pixels from the boundaries of the reference objects with which they were associated. Target objects not colocalized with reference objects were placed at random coordinates within the eroded cell mask, excluding the pixels belonging to reference objects. For this, the reference objects were dilated using a square element with side length = 10 pixels to exclude the colocalization area around each reference object when placing not-colocalized target objects.

(ii) Punctate target and reference. These simulations followed the strategy described in Point i. The only difference was that they did not use the experimental paxillin patch segmentations for reference objects. Instead, reference objects were placed at random coordinates within the eroded cell mask, and a binary image was created where pixels of value one represented the positions of reference objects. These masks were the starting point for the simulations, as described for the previous case.

(iii) Non-punctate target and reference. These simulations also followed the strategy described in Point i. The only difference was that non-punctate target objects were generated by dilating their center points, which were placed near or not near reference objects as in Point i, first with a line of length 7 pixels and a random orientation between 0 and 180° and then with a disk of radius 5 pixels. This generated ellipsoidal target objects with an area of 123 pixels² and eccentricity of 0.79, which were within the range of object size and eccentricity observed for paxillin patches (Fig. S3, A–C).

(i) Punctate target and reference and non-punctate condition. As with the case above of non-punctate reference objects, the paxillin patch segmentations belonging to the employed cell masks were used for the simulations, but this time as condition objects. Based on the input p(TwC) and p(RwC), target and reference objects, were divided into condition-positive and condition-negative groups. Condition-positive reference objects were placed 0.2 pixels away from the boundaries of the condition objects with which they were associated. This very close placement was motivated by practical considerations, making it easier to place target objects (as described next) close to both a reference object and a condition object when needed, as being placed close to one almost guaranteed being placed close to the other. Condition-negative reference objects were placed at random coordinates within the eroded cell mask, excluding the pixels belonging to condition objects after dilating them using a square element with side length = 10 pixels.

To place target objects, they were subdivided further. First, based on p(TwR|TwC), the condition-positive target objects were divided into colocalized or not colocalized with reference objects. Condition-positive target objects colocalized with reference objects were placed 0–2 pixels from condition-positive reference objects. Given the placement scheme for condition-positive reference objects, this almost guaranteed that the target objects were colocalized with condition objects as well. Condition-positive target objects not colocalized with reference objects were placed 0–2 pixels from the boundaries of condition objects, after eliminating the boundary pixels associated with reference objects. Second, based on p(TwR|TnC), condition-negative target objects were divided into colocalized or not with reference objects. Condition-negative target objects colocalized with reference objects were placed 0–2 pixels from condition-negative reference objects, which largely placed them far from condition objects as well. Condition-negative target objects not colocalized with reference objects were placed at random coordinates within the eroded cell mask, excluding the pixels belonging to condition objects and reference objects, both after dilation using a square element with side length = 10 pixels.

(ii) Punctate target, reference, and condition. These simulations followed the strategy described in Point i. The only difference was that they did not use the experimental paxillin patch segmentations for condition objects. Instead, condition objects were placed at random coordinates within the eroded cell mask and a binary image was created where pixels of value one represented the positions of condition objects. These masks were the starting point for the simulations as described in the previous case.

(iii) Non-punctate target, reference, and condition. These simulations also followed the strategy described in Point i. The only difference was that non-punctate target and reference objects were generated by dilating their center points, which were placed near or not near condition objects and/or each other as in Point i, first with a line of length 7 pixels and a random orientation between 0 and 180° and then with a disk of radius 3 pixels for target objects or of radius 5 pixels for reference objects. This generated ellipsoidal target objects with an area of 55 pixels² and an eccentricity of 0.79 and ellipsoidal reference objects with an area of 123 pixels² and eccentricity of 0.79, which were within the range of object size and eccentricity observed for paxillin patches (Fig. S3, A and B).

As mentioned in Results, the conditional colocalization measures employed in this work are conceptually related to the joint probability distribution describing the colocalization of target, reference, and condition. They explore different facets of the joint probability distribution in order to directly address the question of how much the condition influences the extent of target colocalization with reference.

Here we describe a mapping between our conditional colocalization measures and the joint probability distribution. Specifically, suppose that the joint probability distribution for the colocalization of target, reference and condition is given by q(TRC). In this formulation, an object is a “compound object” defined across all three channels, where T = 0/1, R = 0/1, and C = 0/1, indicating that the compound object contains/does not contain a target, reference, and a condition object, respectively. Note that, in order to define compound objects that contain a maximum of one of each object type (T, R, and C), the radius employed for defining compound objects (akin to the colocalization radius employed in our analysis) must be smaller than the minimum inter-object distance within each channel. This puts an upper limit on the radius, which may limit the ability to construct the joint probability distribution in practice, depending on the nature of the analyzed images. Nevertheless, on a conceptual level, there is a mapping between our conditional colocalization measures and the joint probability distribution.

In our simulations, p(TwC), p(RwC), p(TwR|TwC), and p(TwR|TnC) were explicitly defined as input parameters. For the validation of our conditional colocalization analysis, their calculated values were compared directly to their input values. The remaining two conditional colocalization measures, p^rs(Tw(RwC)) and p^rs(Tw(RnC)), and the overall colocalization probability p(TwR), were not explicitly defined. Thus, for validation, their calculated values were compared to their expected values, which were derived from the input p(TwC), p(RwC), p(TwR|TwC) and p(TwR|TnC) and Eqs. 12, 13, 14, 15, 16, 17, 18, 19, 20, 21, and 22 above.

Eq. 23 reflects the law of total probability relating p(TwR) to the conditional probabilities p(TwR|TwC) and p(TwR|TnC).

To reduce image SNR from its original value to 4, 3, and 2 (Fig. 4), Gaussian white noise was added to the images. The standard deviation of the added noise was adjusted to reach the desired SNR as follows:

If the introduction of additional noise rendered some intensity values negative, the mean background intensity was increased to ensure that all intensity values remained positive after the added noise.

Fig. S1 shows example illustrating the randomization procedure for non-punctate objects and convergence of calculated coincidental measures upon randomization. Fig. S2 Illustrates the colocalization distance for two general, non-punctate objects and its dependence on their relative location, size and shape. Fig. S3 shows properties of segmented paxillin patches and further colocalization analysis of simulated non-punctate objects. Fig. S4 shows the results of applying conditional colocalization analysis to an experimental negative control. Fig. S5 shows the results of applying conditional colocalization analysis to non-punctate images of activated β₁-integrin, pFAK and paxillin.

We thank Dr. Tieqiao Zhang for microscopy support, Dr. Anthony R. Vega for initial programming of code for pairwise colocalization of punctate objects, and Dr. Nicolas Touret (University of Alberta, Edmonton) for critical reading of the manuscript.

This work was supported by funding from National Institutes of Health/National Institute of General Medical Sciences (R35 GM119619), the Welch Foundation (I-1901), the National Science Foundation (MCB-2114417), and the University of Texas Southwestern Endowed Scholars Program to K. Jaqaman. J. Vega-Lugo was an honorary trainee of the National Institutes of Health Molecular Biophysics Training Grant (5T32GM131963; PI: Dr. Yuh Min Chook) and was also supported by a Diversity Supplement to National Institutes of Health grant R35 GM119619 (PI: K. Jaqaman).

The authors declare no competing financial interests.

Author contributions: K. Jaqaman, J. Vega-Lugo and B. Rocha-Azevedo designed research; J. Vega-Lugo and K. Jaqaman developed algorithm; J. Vega-Lugo wrote software; B. Rocha-Azevedo and A. Dasgupta performed experiments; J. Vega-Lugo and B. Rocha-Azevedo performed analysis; K. Jaqaman, J. Vega-Lugo and B. Rocha-Azevedo wrote paper with input from A. Dasgupta.