KMT5C encodes robust heterochromatin retention and liquid-like behavior using limited sequence features

Cells use multiple strategies to compartmentalize functions through a combination of membrane- bound and membraneless organelles. The latter represent complex assemblies of biomolecules that coalesce into a dense phase through low affinity, multivalent interactions and undergo rapid exchange with the surrounding dilute phase. We describe a liquid-like state for the lysine methyltransferase KMT5C characterized by diffusion within heterochromatin condensates but lacking appreciable nucleoplasmic exchange. Retention was strongly correlated with reduction of condensate surface area, suggesting formation of a liquid droplet with high connectivity. This behavior mapped to a discrete domain whose activity was dependent on multiple short linear motifs. Moreover, it was strikingly resilient to marked phylogenetic differences or targeted changes in intrinsic disorder, charge, sequence, and architecture. Collectively, these findings show that a limited number of sequence features can dominantly encode multivalency, localization, and dynamic behavior within heterochromatin condensates to confer protein retention without progression to a gel or solid.


SUMMARY
Cells use multiple strategies to compartmentalize functions through a combination of membranebound and membraneless organelles. The latter represent complex assemblies of biomolecules that coalesce into a dense phase through low affinity, multivalent interactions and undergo rapid exchange with the surrounding dilute phase. We describe a liquid-like state for the lysine methyltransferase KMT5C characterized by diffusion within heterochromatin condensates but lacking appreciable nucleoplasmic exchange. Retention was strongly correlated with reduction of condensate surface area, suggesting formation of a liquid droplet with high connectivity. This behavior mapped to a discrete domain whose activity was dependent on multiple short linear motifs. Moreover, it was strikingly resilient to marked phylogenetic differences or targeted changes in intrinsic disorder, charge, sequence, and architecture. Collectively, these findings show that a limited number of sequence features can dominantly encode multivalency, localization, and dynamic behavior within heterochromatin condensates to confer protein retention without progression to a gel or solid.

INTRODUCTION
Cells contain membraneless organelles of varying complexity where constituents partition between a dense liquid demixed phase and a dilute surrounding phase (1, 2). Formation of these biological condensates reflects low affinity, multivalent interactions that allow proteins to coalesce into distinct aqueous compartments (2). Despite the stable appearance of such assemblies, their constituents rapidly exchange with their surroundings (3,4). Nevertheless, they can adopt a range of dynamic and material states, including gels or solids (5)(6)(7). This can be exemplified by the nucleolus, which is made up of three biophysically distinct phase separated compartments (8), each with specialized functions (9)(10)(11). At present, however, we have an incomplete understanding of the range of dynamic behaviors exhibited in the context of biomolecular condensates and how this relates to their physiological roles.
Condensate formation can be conceptualized by categorizing components as either scaffolds or clients (12)(13)(14). Scaffold molecules are essential to initiating phase separation and provide a platform for recruitment of client proteins. The valency or number of interaction sites on scaffold and client molecules defines network connectivity and the material properties of the assembly (15,16). In this regard, proteins with higher valency will phase separate at lower concentration and form condensates with greater stability (15,16). Although the sequence features that sustain these interactions have been defined in many cases, they display considerable diversity (2,(17)(18)(19). These include complex interaction types mediated by single amino acids, short linear motifs, regions of intrinsic disorder, structured domains, and nucleic acids (17,20,21). Given the repertoire of proteins that partition to condensates is not complete, there are opportunities to discover new scaffold or client proteins with unique sequence and biophysical properties, as well as those that use potentially novel strategies.
In the present study, we assessed inter and intra-chromocenter mobility for KMT5C in comparison to CBX5 (28, 29) and MeCP2 (42)(43)(44), which provided known benchmarks for heterochromatin phase separation. Surprisingly, we found that KMT5C sampled the entire chromocenter volume on the time scale of seconds but underwent only minimal exchange with the surrounding nucleoplasm over the same period. In contrast, both CBX5 and MeCP2 exchanged rapidly, as previously reported (36, [45][46][47]. The unique dynamic behavior of KMT5C was associated with robust chromocenter partitioning and reduction of surface area. Together, this suggests the formation of a liquid-like condensate that effectively sequestered KMT5C within chromocenters without imposing the same constraint on CBX5 or MeCP2. This activity was consolidated within a discrete KMT5C protein segment that displayed considerable sequence divergence amongst tetrapods. Nevertheless, its localization and dynamic behavior were largely conserved and conferred by a unique combination of short linear motifs that functioned in a variety of sequence contexts. Overall, our work provides a model to understand how sparse sequence features can confer efficient protein retention and liquid-like behavior within a biological condensate.

KMT5C defines a novel dynamic state in phase-separated constitutive heterochromatin
In fluorescence recovery after photobleaching (FRAP), CBX5 rapidly exchanged between chromocenters and the nucleoplasm (36, 37, 46,47), whereas SUV39H2 and KMT5C appeared to be immobile when entire chromocenters were bleached (35-37). Although concluding the latter two proteins created a stable scaffold, these studies did not assess mobility within chromocenters, which is required to determine if they adopt a liquid or solid-like state. For our analyses, we used the immortalized mouse mammary epithelial cell line NMuMG (48), which showed prominent H4K20me3 staining that was largely confined to chromocenters (fig. S1A), indicating the underlying KMT5C methylation pathway was functional. Importantly, this limited H4K20me3 distribution was also observed in epithelial cells of the mouse pubertal mammary gland ( fig. S1B), establishing NMuMG cells as a relevant model system to understand how KMT5C dynamics relate to the distribution of its enzymatic product.
To evaluate the possibility of KMT5C mobility within constitutive heterochromatin, we queried KMT5C dynamics in mouse chromocenters in comparison to CBX5 and MeCP2, which exhibit pericentromeric enrichment, nucleoplasmic exchange (36, 45), and undergo liquid-liquid phase separation (28,29,(42)(43)(44). As previously shown with full chromocenter bleaching (35-37), KMT5C displayed long recovery times, indicating minimal exchange with the surrounding nucleoplasm ( Fig. 1A, Total). Upon partial bleaching, however, fluorescence recovered (movie S1) and progressed from the non-bleached portion of the chromocenter (Fig. 1A, Partial). This occurred on a timescale where there was no appreciable recovery of fully bleached chromocenters in the same nucleus (Fig. 1A, Total), establishing that KMT5C moved readily within chromocenters but did not efficiently exchange. In other words, KMT5C equilibrated within individual chromocenters with minimal influx from the surrounding nucleoplasm. Both MeCP2 and CBX5, however, showed near complete recovery in fully bleached chromocenters over the same time, which was readily apparent in both imaging (Fig. 1B) and FRAP analysis ( Fig. 1C and movies S2 and S3). KMT5C therefore appeared to be effectively confined to individual chromocenters when compared to MeCP2 or CBX5 (Fig. 1C, schematics) and displayed significantly higher partitioning to chromocenters than either protein (Fig. 1D). Thus, KMT5C, MeCP2, and CBX5 define distinct dynamics states in constitutive heterochromatin defined by progressively increasing intra-chromocenter mobility and nucleoplasmic exchange that is inversely correlated to partitioning (Fig. 1E).
The dynamic behavior of KMT5C was reinforced using two additional imaging approaches. First, photobleaching half of the nucleus established the pan-nuclear nature of  fig. S2C, inset). In the context of phase separation, the behavior of KMT5C signified a remarkable barrier to exit that has not been observed for other proteins without the transition to a solid or gel state, which would render them immobile (49,50).

KMT5C exhibits liquid-like properties in constitutive heterochromatin condensates
The dynamic behavior of KMT5C suggested retention within a liquid droplet. Spot bleaching of a defined volume within individual chromocenters ( Fig. 2A) provided clear evidence of internal mobility with a t1/2 of ~4.4s and an apparent diffusion (Dapp) constant of 1.9 x 10 -2 µm 2 s -1 (Fig.   2B) (51). Internal mixing was also evident in convergence curves (52)  Bisection of chromocenters with laser micro-irradiation led to rapid KMT5C retraction from the damaged area to create two lobes that preserved the retentive behavior of the original domain by FRAP analysis (Fig. 2C, fig. S3C, and movies S4 and S5), indicating a lack of exchange between lobes. The rounded shape of these compartments provided further evidence of liquidity where KMT5C appeared to drive minimization of surface area (1). This was not apparent with either CBX5 or MeCP2, which were characterized by distinct responses. Although CBX5 initially showed depletion from the damage zone, this was transient and accompanied by redistribution throughout the damaged and non-damaged compartments, where it continued to show rapid recovery after photobleaching (Fig. 2C, fig. S3C, and movies S6 and S7). The localization and mobility of MeCP2, however, appeared unperturbed by laser micro-irradiation (Fig. 2C, fig.  S3C, and movies S8 and S9). Lastly, KMT5C was associated with a significantly higher degree of sphericity compared to either MeCP2 or CBX5 in transfected cells (Fig. 2D, p ≤ 0.001), indicating it was more effective at reducing surface area at the boundary between the chromocenter and nucleoplasm. KMT5C therefore demonstrates that phase separation can achieve robust protein sequestration while nevertheless retaining liquid-like internal mobility in the context of heterochromatin-associated condensates.

A rudimentary heterochromatin targeting domain underlies CRD activity
At the primary sequence level, KMT5C comprises an enzymatic SET domain that catalyzes Inclusion of the adjacent RBD, which confers nucleolar localization during quiescence (54), did not alter localization or dynamics (not shown). It therefore appears that the interactions underlying the dynamic behavior of full-length KMT5C are largely contained within this discrete region, which we will refer to as the Chromocenter Retention Domain (CRD) to reflect its activity.
We sought to understand the overall architecture and driver(s) of CRD behavior through examination of sequence features. Nevertheless, sequence identity amongst mammalian CRDs was limited to 19 of 57 residues and there was considerable variation within individual ordersfor example, 23 positions were variable within Primates. Moreover, manual curation was required to identify non-mammalian CRD counterparts, which were limited to bird and reptile KMT5C orthologs. Notwithstanding marked divergence, we identified three conserved features: a His-x4-Cys-x3-4-Cys motif and a L-x5-6-V-x1-L hydrophobic motif that were always separated by a single amino acid (Fig. 4A). These were followed by a related L-x5-V-x1-L motif that was embedded within distinct blocks of homology in each taxonomic class and separated by segments of highly variable length and sequence composition ( Fig. 4A and table S1). Homology within the CRD therefore occurred as two discrete segments, designated as conserved domains 1 and 2 (CD1 and CD2), separated by a non-conserved linker region, suggesting a bipartite structure. This fit with previous analyses of the mouse KMT5C region spanning residues 350-412 where each half could independently confer chromocenter localization, albeit with lower efficiency and elevated mobility in FRAP experiments (35). We confirmed this finding with full and partial chromocenter bleaching of CD1 and CD2-mEmerald fusion proteins, which identified CD1 as the main driver of CRD behavior ( fig. S4A and S4B).
To query if differences in sequence, physicochemical properties, or architecture affect the localization and dynamic behavior of the CRD, we assessed versions from Homo sapiens (Hs),  (Fig. 4B). The fact that these behaviors were shared by all orthologs indicates they are driven by a limited number of sequence features that are modulated by species-specific differences.

CRD localization and dynamics are dominantly encoded by Short Linear Motifs (SLiMs)
We examined the role of 4 highly conserved motifs in the CRD, focusing on how C 362 C 366 , Mutation of individual motifs (Fig. 5A) did not abrogate chromocenter localization but caused varying degrees of nucleoplasmic distribution (Fig. 5B). All mutants were characterized by elevated mobility in both full and partial bleaching conditions ( Fig. 5C and movies S11-S15).
This was most pronounced for C 362 C 366 , H 357 H 365 , and L 368 V 374 L 376 where the partial and full recovery curves were largely superimposable (Fig. 5C), indicating differences between internal and external mobility were abrogated (49). In contrast, L 398 V 404 L 406 had a modest effect that was characterized by similar increases in internal mobility and exchange. The nature of differences between L 368 V 374 L 376 and L 398 V 404 L 406 further support CD1 making a greater contribution to CRD localization and dynamics ( fig. S4A and S4B). Lastly, we assessed the mammalian-specific motif, W 359 W 390 W 392 Y 396 , because it resembled the aromatic cage motifs found in chromatin reader modules that bind histone methyl-lysine residues (55)(56)(57). This mutant also ameliorated the difference between partial and full chromocenter fluorescence recovery but had an intermediate effect when compared to C 362 C 366 , H 357 H 365 , and L 368 V 374 L 376 (Fig. 5C). Each motif therefore constituted a distinct contribution to CRD localization and dynamics, and no single motif was absolutely required for these activities. Together, the results indicate that multiple Short Linear Motifs (SLiMs) act in a highly synergistic manner to confer KMT5C chromocenter retention.
Intrinsically disordered regions (IDRs) play an important role in protein liquid-liquid phase separation (58-60) but were highly variable for the CRD between and within taxonomic classes ( Fig. 6A and table S2). Specifically, while each CRD showed a low propensity for disorder in CD1 that was minimal over the H-C-C segment (Fig. 4A), there was considerable variability in disorder profiles for the linker and CD2 regions (Fig. 6A, top). In fact, the only taxonomic class that displayed a conserved disorder profile was birds, where CD1 and CD2 were separated by a highly disordered linker region regardless of length (Fig. 6A, bottom, and table S2). To assess the potential contributions of intrinsic disorder to the dynamic properties in a controlled manner, we used an IDR series that is known to modulate phase separation behavior when linking globular domains (61). These 40 amino acid segments were used to replace the 6 amino acid linker in the mouse CRD to create chimeras with an extended region of disorder ( Fig.   6B and table S2) that varied 8-fold in the fraction of charged residues, while net charge per residue varied from -0.213 to 0.085 ( Fig. 6C and table S1). Although each CRD-IDR chimera conferred efficient chromocenter localization ( fig. S5), all five showed a higher level of fluorescence recovery than the wild-type CRD upon partial bleaching (Fig. 6D). At the same time, the CRD-IDR1 to IDR4 proteins produced recovery curves that were superimposable with the wild-type CRD upon full chromocenter bleaching (Fig. 6D). Only CRD-IDR1, which had the highest charge fraction and net charge per residue (Fig. 6C), showed a modest increase in exchange (Fig. 6D). Importantly, these CRD-IDR chimeras establish that exchange and internal mobility are separable. Moreover, together with the mutational analyses (Fig. 5), they indicate that SLiMs dominantly encode the retentive behavior of the CRD, while differences in disorder between CRDs may tune material properties within chromocenter condensates.

DISCUSSION
Several studies have assessed KMT5C/SUV420H2 mobility in the context of mouse chromocenters, with all concluding it was immobile (35-37). Importantly, this conclusion was also reached using a knockin allele that expressed a KMT5C-GFP fusion protein at endogenous levels (35), indicating it was not a consequence of exogenous expression. Revisiting dynamics in the context of heterochromatin condensate formation (24, 62), we found that KMT5C and the CRD displayed liquid-like behavior within chromocenters, but underwent limited exchange with the surrounding nucleoplasm (Figs. 1-3). This suggests that the dynamic behavior of KMT5C and the CRD is shaped by strong chromocenter avidity, which could reflect differences in solvation properties between the dilute and dense phases (63,64) that drive robust protein compartmentalization without gelation. This behavior was fundamentally different from that of MeCP2 and CBX5 (Fig. 1), indicating phase separated heterochromatin domains can accommodate a continuum of dynamics states.
Liquid-liquid phase separation has been increasingly evoked to account for the numerous molecular assemblies within the nucleus, involving a diverse range of size scales and regulatory processes (22,65,66). Whether all are bona fide or can be described by other biophysical processes remains to be determined (4,17,65,67). This applies to the heterochromatin assemblies under study here, which have been described as liquid-liquid phase separated (25, 28, 29, 62) and alternatively as "collapsed polymer globules" (24, 52, 68). In the latter case, the dynamic behavior of proteins would be dictated by percolation through the chromatin polymer (52). The preferential internal mobility of KMT5C within chromocenters and reduction of surface area (Figs. 1 and 2), however, are consistent with liquid-liquid phase separation (52).
Another key consideration in interpreting these dynamic behaviors is the material state of the underlying chromatin scaffold in terms of being liquid-like or solid-like (26, 27). Notably, KMT5C exhibited internal remixing within chromocenter condensates even when the underlying chromatin did not (26)-this solid-like chromatin behavior was also seen in studies of CBX5 and MeCP2 chromocenter mobility (52). As a result, KMT5C movement within chromocenters and the absence of nucleoplasmic exchange do not reflect stable binding (35-37) to mobile, liquidlike chromatin (27). Lastly, the positional stability conferred by this chromatin material state likely limits the opportunity for chromocenters to undergo fusion and/or fission, which is observed for condensates that form by liquid-liquid phase separation (17,69,70). Nevertheless, insomuch as laser micro-irradiation induced a phenomenon akin to fission, KMT5C exhibited the expected behavior for a liquid-like droplet (Fig. 2).
Proteins that undergo phase separation can be classified as either scaffolds that initiate condensate formation or clients that modulate condensate properties through low-affinity interactions (12,15,69). CBX5 (71) and MeCP2 (72)  For a given scaffold, condensate stability can be increased by client proteins that have higher valency or those that add to connectivity without disrupting scaffold-scaffold interactions, both of which support elevated partitioning of the client protein to the condensate (15). In contrast, low valency clients reduce condensate stability (15) by using capping mechanisms that attenuate or block network connections and thereby oppose phase separation (74). Of these two scenarios, the robust partitioning of KMT5C and the CRD are consistent with creating a high level of connectivity. Likewise, that they encompassed dynamic behaviors resembling viscous liquids with respect to internal mobility, or viscoelastic materials with respect to exchange between phases in FRAP analyses (50), requires interactions of sufficient number and strength to modulate emergent properties (12,15). The fact that these features were apparent upon exogenous over-expression indicates connectivity persists over a wide concentration range (15) and that the underlying scaffold is not limiting or that KMT5C can undergo client to scaffold switching.
The CRD provides a minimalistic model to understand sequence features that support condensate retention within an otherwise liquid-like state. Whereas CBX5 (28, 29) and MeCP2 (42,43,75) use a combination of structured domains (stickers) and intrinsically disorder regions (spacers) that contribute to specific and non-specific interactions (2,61), the dynamic and partitioning behaviors of KMT5C have been consolidated into the CRD. In comparison to CBX5 and MeCP2, which are highly conserved, it was also clear that the CRD was under very different evolutionary constraints. Comparing this domain across tetrapods revealed limited conservation that was only apparent after manual alignment, but nevertheless conferred efficient chromocenter localization ( fig. S4C). It was not until the introduction of discrete mutations in a series of conserved motifs that the CRD showed a gradual shift towards full recovery through nucleoplasmic exchange, indicative of progressive weakening or loss of interactions within the condensate (15). Strikingly, the CRD imparted these behaviors in a broad range of sequence contexts and in the presence of multiple heterologous protein sequences that substantially varied physicochemical properties and disorder content. CRD phase separation therefore appears to be dominantly encoded by a collection of SLiMs that can act as strong attraction promoters (76) to support robust retention, while sequence differences tune dynamic properties within chromocenters (Fig. 6) (77). In the context of protein features or molecular grammar that drive phase separation (18-20, 59, 78), the CRD appears to be an outlier over and above its unique dynamic properties.
Pericentromeric heterochromatin is the principal target for KMT5C-mediated catalysis of H4K20me3 (33, 34, 38, 39) (fig. S1). As a result, the retention of KMT5C within mouse chromocenters coincides with its product distribution and provides a means for the strict spatial control of this epigenetic mark. Its requirement for pre-existing H3K9me3 and HP1 chromobox proteins (33, 37) resembles the dependence of H2A lysine-119 ubiquitylation by RNF2 on formation of a phase separated compartment driven by binding of the CBX2 chromobox protein to H3K27me3 (79). Here, condensate formation is thought to provide a "reaction crucible" for localized activity of RNF2 (79). Consistent with this, perturbation of KMT5C chromocenter localization in Suv39h1/2 knockout cells is associated with spreading of H4K20me3 outside of constitutive heterochromatin (33). Retention also provides a means to limit KMT5C access to substrates, which can account for why the trimethyl modification represents only a small fraction of H4K20 methylation states in mouse (35) and human (40, 41) cell lines. This differs from KMT5B catalysis of H4K20me2 (34), which accounts for the majority of H4K20 methylation states in mouse (35) and human (40, 41) cell lines and is therefore broadly distributed across the genome. At the DNA level, this reflects the association of H4K20me3 with arrays of tandem satellite repeats that make up the bulk of constitutive heterochromatin at pericentromeres (39, [80][81][82]. Although these sequences are highly divergent between species (81,(83)(84)(85), they are thought to share a core set of epigenetic features (86,87). Our findings with divergent CRD orthologs support a model wherein they recognize constitutive heterochromatin features that are conserved across distant evolutionary timescales. At the same time, the CRD may accommodate differences in rapidly evolving pericentromeric repeat sequences (84,88) and it will be interesting to determine if and how orthologs contribute to differences in constitutive heterochromatin between and within species (30, 89-91).

Cell Culture and transfection
Cells were cultured at 37°C and 5% CO2 in a humidified incubator. All cell lines were grown in DMEM containing 10% FBS. D5 (Suv39h1/2 knockout) and W8 (Suv39h1/2 wild-type) mouse embryonic fibroblast cells lines (32) were obtained from Dr. Thomas Jenuwein. All other analyses were carried out using the mouse NMuMG breast cell line. Cells were transfected by lipofection using Effectene (Qiagen) 1 day prior to experiments. Expression plasmids were synthesized (www.biomatik.com), obtained from the addgene repository (www.addgene.org), or previously described (92). Pubertal mouse mammary gland sections were provided by Dr. Ing Swie Goping in accordance with the guidelines and regulations set forth by the Canadian Council on Animal Care and approved by the University of Alberta Health Sciences 2 Animal Care and Use Committee (Protocol #AUP00000386). Chromocenter spot-bleaching was performed on a PerkinElmer UltraVIEW ® spinning disc confocal microscope with Volocity ® 6.3 software. Data was acquired using a 100x 1.4 NA Zeiss Plan-APOCHROMAT objective and a 1.6x optovar at 37°C with humidification and 5% CO2. The center of each chromocenter was selected using the Volocity ® Spot Tool and bleached using 120 cycles of the 488 nm laser at maximum power. Image acquisition was set to record 5 pre-bleach frames followed by 40 post-bleach frames at 1-second intervals. The recorded time series were imported to ImageJ for registration. Mean intensities were acquired for the bleached chromocenter, a reference chromocenter, and background. Individual curves were normalized in Microsoft Excel, averaged, and plotted with GraphPad Prism using a mono-exponential model (f(t) = A(1-e (-t/t) ) for curve fitting. The radius of bleach spots was determined using single images or stacks of paraformaldehyde-fixed KMT5C-mEmerald expressing cells by measuring the full width at half maximum intensity and dividing the average by two. The apparent diffusion constant (0.019 µm 2 s -1 ) was calculated from the measured radius (~0.32 ±0.08 µm) of the bleached spot and the time value (t) that indicated full fluorescence recovery using the following formula: Dapp = r 2 /t (51).

Immunofluorescence
Cells grown on adherent coverslips were fixed in 4% paraformaldehyde in PBS for 10 min,

Image analysis
Images were deconvolved with Huygens Professional version 19.04 (Scientific Volume Imaging, http://svi.nl). Z-stacks were imported into Imaris 9.3 (Oxford Instruments) and cropped to generate 3D images of single nuclei. The Imaris surfaces function was used to generate masks based on intensity to encompass all chromocenters and the statistics function was used to measure the volume, sphericity, and total number of chromocenters for output to a Microsoft Excel spreadsheet. For partition coefficients, line scans through nuclei of undeconvolved 3D images were recorded using ImageJ. For each cell, the partition coefficient was calculated by subtracting the background level from the maxima of the brightest chromocenter and dividing it by the background corrected fluorescence intensity of the nucleoplasm. For each protein analyzed, measurements were taken from n>10 nuclei. Kymographs were acquired for regions of interest from laser micro-irradiation time series using the ImageJ Multi Kymograph tool.

Sequence analyses
In addition to recovery of all tetrapod KMT5C orthologs annotated in the NCBI database (https://www.ncbi.nlm.nih.gov/gene/84787/ortholog/?scope=32523&term=KMT5C), the Basic Local Alignment Search Tool (BLAST) (94) was used to expand this list to include unannotated orthologs in birds and reptiles using tblastn and a range of genome and transcript level sequence databases. Multiple alignments were carried out independently for CD1 and CD2 using ClustalW (95). Disorder profiles were derived using PONDR â (http://www.pondr.com/) and VL-XT mode.

Statistical analysis
For partition coefficients, significance was evaluated using the Kruskal-Wallis one-way analysis test and individual comparisons between proteins were done using the Wilcoxon rank sum test.
Statistical significance of sphericity data was evaluated using the embedded stats function within PlotsOfDifferences, which calculates p-values using a randomization test (99).

21.
H. X. Zhou, V.      Partial (open circles) conditions. The wild-type CRD (from figure 4) is shown for reference.