<figure><img src="image_1.png" /></figure>

HAL
open science

# Scaling characteristics of modelled tropical oceanic rain clusters

Chee-kiat Teo, Tieh-yong Koh, Kevin Cheung, Bernard Legras, Hoai-nguyen Huynh, Lock-yue Chew, Leslie Norford

► To cite this version:

Chee-kiat Teo, Tieh-yong Koh, Kevin Cheung, Bernard Legras, Hoai-nguyen Huynh, et al. Scaling characteristics of modelled tropical oceanic rain clusters. Quarterly Journal of the Royal Meteorological Society, 2021, 10.1002/qj.3959. hal-03142906

HAL Id: hal-03142906

https://hal.science/hal-03142906v1

Submitted on 23 Oct 2025

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_2.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
RESEARCH ARTICLE

# Scaling characteristics of modelled tropical oceanic rain clusters

Chee-Kiat Teo<fnref n="1" />

Tieh-Yong Koh<fnref n="2" />

Kevin K. W. Cheung<fnref n="3" />

<figure><img src="image_3.png" /></figure>

Bernard Legras<fnref n="4" />

Hoai-Nguyen Huynh<fnref n="5" />,<fnref n="6" />

Lock-Yue Chew<fnref n="7" />

Leslie Norford<fnref n="8" />

**Correspondence**

K. Cheung, Climate and Atmospheric
Science, NSW Department of Planning
Industry and Environment, Parramatta,
NSW, Australia.

Email: kkwc09@gmail.com

**Funding information**

National Research Foundation Singapore

**Abstract**

The scaling exponents of the distributions of cluster rain amount, *R*, and cluster size, *A*, for oceanic rain clusters over the Indian and Pacific warm pools, and the intertropical convergence zones over the eastern Pacific and the tropical Atlantic, were obtained from a set of regional climate model downscaling products. The main aim of the investigation is to compare the model cluster's scaling characteristics with those obtained from observations that have been reported previously. The scaling exponents for the model were found to be different across the ocean basins indicating the lack of universality in the modelled rain cluster distributions. The scaling exponent for the conditional mean of *R* given *A* = *a*, E(*R*|*a*), was found to be the same across the different ocean basins, and the estimated value of the exponent agrees with that obtained from satellite-observed rain clusters. However, no crossover in the scaling of E(*R*|*a*) in the model for cluster size larger than mesoscale was seen, unlike those reported elsewhere based on observations. The implication is that in the model the intensification of rain with cluster size continues up to synoptic scale. Through simple scaling arguments it is believed that before the break in scaling that has been identified from an observational study elsewhere, the model simulates the fundamental mesoscale dynamics well and thus estimated the E(*R*|*a*) in agreement with observations. Whether such a difference in transition of scaling between the modelled and observed rain cluster scaling behaviour depends on the model details, for example the convection parametrization, needs further clarification.

**KEYWORDS**

1. tools and methods: Regional and mesoscale modelling, 2. scale: global, mesoscale, 3. physical phenomenon: clouds, convection, rainfall

<footnote n="1">¹Centre for Climate Research Singapore, Meteorological Services Singapore, Singapore</footnote>

<footnote n="2">²Centre for University Core, College of Lifelong and Experiential Learning, Singapore University of Social Sciences, Singapore</footnote>

<footnote n="3">³Climate and Atmospheric Science, Science Economics and Insights Division, NSW Department of Planning Industry and Environment, Parramatta, New South Wales Australia</footnote>

<footnote n="4">⁴Laboratoire de Météorologie Dynamique, CNRS/ENS-PSL/Sorbonne U./Ecole Polytechnique, Paris, France</footnote>

<footnote n="5">⁵Institute of High Performance Computing, Agency for Science, Technology and Research, Singapore</footnote>

<footnote n="6">⁶Complexity Institute, Nanyang Technological University, Singapore</footnote>

<footnote n="7">⁷School of Physical and Mathematical Sciences, Nanyang Technological University, Singapore</footnote>

<footnote n="8">⁸CENSAM, Singapore-MIT Alliance for Research and Technology, Singapore</footnote>
# 1 | INTRODUCTION

Several investigations have demonstrated that deep convective rain possesses universal behaviour across the different tropical ocean basins and as a result, the probability distributions of the rain clusters, when suitably normalised, are identical (Peters and Neelin, 2006; Peters et al., 2009; Teo et al., 2017). One approach to explain how such universality can arise is to appeal to the hypothesis that tropical deep convective rain is a self-organised criticality (SOC) among the other models (Traxl et al., 2016). Such a view was first forwarded by Andrade et al. (1998) and then more formally proposed by Peters and Neelin (2006) where they pointed out that the tropical convective atmosphere has the following ingredients for SOC: (a) the existence of a column water vapour content threshold beyond which the convective rain intensity increases exponentially, and (b) the clear separation of the relaxation (convective) time-scale and the driving time-scale of the (relatively slower) build-up of the atmospheric water vapour content from the latent heat fluxes from the ocean surface and large-scale moisture advection.

It has been proposed that the ubiquitous presence of deep convective rain clusters of a wide range of spatial scales over the tropical oceans is yet another aspect of the tropical convective rain as an SOC (Peters and Neelin, 2006; Neelin et al., 2008; Teo et al., 2017, abbreviated as T17 hereafter). T17 suggested that the rain clusters are analogous to avalanching sandpiles – the archetypal SOC model (e.g. Bak et al., 1987), with the rain cluster size, A, and the cluster rain amount, R (total rain intensity over the cluster area), corresponding to the avalanche area and avalanche size of a critical sandpile. The probability density distributions for both A and R are by virtue of being an SOC assumed to follow simple scaling (Christensen and Moloney, 2005, p. 273):

$$f_S(s) \sim s^{-\zeta_S} \mathcal{G}_S \left( \frac{s}{s_c} \right), \quad s > s_0; s \in \{A, R\} \quad (1)$$

where $f_S(s)$ is the probability density distribution, $\zeta_S$ is the scaling exponent, and $\mathcal{G}_S$ is the scaling function expressing the departure from a strict power law due to the finite size of the system (the oceanic domain where these tropical mesoscale rain clusters organise; see Figure 1). This scaling function $\mathcal{G}_S(y)$ has the following property: it decays faster than a power law for $y \gg 1$, and stays (almost) as a constant for $y \ll 1$. Hence the probability distribution exhibits a strict power-law regime ($f_S(s) \sim s^{-\zeta_S}$) only for values of s larger than the lower cut-off, $s_0$, but much smaller than the upper cut-off, $s_c$. For simple scaling to be valid, the upper cut-off, $s_c$, increases monotonically with the system size, while the lower cut-off, $s_0$, is independent of the system

size. By analysing the distributions of satellite-observed rain clusters over oceanic tropical rain clusters, T17 concluded that the tropical oceanic rain clusters are universal (i.e. consistent across the ocean basins) with the scaling exponents $\zeta_A$ and $\zeta_R$ being 5/3 and 3/2 respectively. T17 also found two scaling regimes for conditional mean of R given $A = a$, $E(R|a)$ from the observed rain clusters: for a rain cluster of length-scale smaller than ~320 km, $E(R|a) \sim a^{4/3}$. For a cluster with length-scale larger than 320 km, $E(R|a)$ was proportional to $a$. They suggested that such scaling behaviour could mean that rain clusters tend to grow or decay nonlinearly below 320 km while larger clusters of size beyond 320 km grow into or decay from superclusters without intensification.

Model evaluation based on scaling of the modelled variables has been emerging (e.g. Skamarock, 2004; Kahn et al., 2011; O'Brien et al., 2013; Fang and Kuo, 2015). In this work, we analysed the scaling characteristics of oceanic tropical rain clusters from a regional climate model. The main questions we wish to address in this work are: (a) whether the distributions of the model rain clusters are universal across all the ocean basins, and (b) whether the model can reproduce the same scaling behaviour as what has been observed for actual rain clusters. The resulting analyses serve as a set of diagnoses for the self-organising behaviour of model rain clusters in the tropics. Increasingly, multiscale interactions are recognised to be an important aspect of the dynamics of tropical weather systems from mesoscale convective complexes and synoptic-scale convectively coupled equatorial waves, up to planetary-scale systems like the Madden-Julian Oscillation (MJO) (e.g. Nakazawa, 1988; Biello and Majda, 2005; Tulich and Mapes, 2008; Kiladis et al., 2009). Therefore, diagnoses that encapsulate the scaling behaviour of rain clusters, such as those presented in this work, are indicative of the model's performance in the tropics.

At the same time, investigating scaling behaviour of rain clusters from the model may be a step towards identifying the essential mechanisms that explains the emergence of SOC in the self-organisation of tropical rain clusters. Besides the model proposed by Stechmann and Neelin (2014) for the temporal rain clusters that exhibited SOC, the spatial extension of a similar model in Hottovy and Stechmann (2015) and the physically based stochastic model in Ahmed and Neelin (2019), the authors are not aware of any other physically based SOC model that accounts for the system size distributions of self-organising deep convection clusters. Obtaining a representative set of scaling exponents of the model rain cluster distributions is therefore a useful first demonstration of similar analyses of rain clusters from different model physics configurations or from various other
<figure><img src="image_4.png" /></figure>

<figure><img src="image_5.png" /><figcaption>**FIGURE 1** (a–d) Selected domains for different Z (domain size, the same in both directions) together with the mean model rain intensity. Boxes with solid (dotted) edges are domains that passed (failed) the selection criteria. The contours are the mean model rain intensity (mm·hr⁻¹) from 1988 to 2015 in logarithmic scale. The two boxes with thicker borders in the AO domain in (b) are two additional Z = 48 domains added for analysis (see main text for details)</figcaption></figure>

models, such as those from Quinn and Neelin (2017a; 2017b). Ideally, some current or future meteorological models or model physics configurations ought to reproduce the observed scaling exponents and shed light on the actual physical-dynamical mechanisms that give birth to SOC behaviour in tropical rainfall. This would facilitate the understanding of SOC in tropical rainfall away from the conventional representations of SOC which are mostly “toy” mathematical representations not based on the universal laws of physics.

The remainder of the article is as follows: Section 2 outlines the model products used and how the rain clusters were obtained. Section 3 describes how the

cluster probability density distributions were estimated and reports on their general characteristics. A more detailed analysis of the scaling exponents of the model rain clusters is in Section 4. Further remarks and discussions on the implications of our findings are in Section 5 before concluding in Section 6.

## 2 | DATA AND PROCESSING

Twenty-eight years of gridded three-hourly model rainfall accumulation from April 1988 to March 2015 produced with the Weather Research and Forecasting model (WRF:
Skamarock *et al.*, 2008) were made available for this investigation. The rainfall is the sum of the rainfall from the cloud microphysics scheme and the parametrized cumulus convection. A brief description of the model is summarized below. Details of the model configurations and model performance can be found in Fonseca *et al.* (2019).

The model domain was a tropical channel confined within 36°S–40°N. The horizontal resolution of the model was 36 km, with 37 vertical levels up to 30 hPa that were spaced more closely together in the planetary boundary layer (PBL) and the tropical tropopause layer. The model was initialized with Climate Forecast System Reanalysis (CFSR: Saha *et al.*, 2010) 6-hourly data. Analysis nudging of several model fields towards the CFSR values was applied in the model runs (Fonseca *et al.*, 2019): the model horizontal winds and potential temperature perturbations were nudged in the lower stratosphere to provide observed upper boundary conditions associated with the Quasi-Biennial Oscillation (Andrews *et al.*, 1987), while the water vapour mixing ratio was nudged in the mid-troposphere to maintain a close correspondence to observed tropical moisture distribution in the model runs, each of uninterrupted one-year duration, for a good simulation of the MJO (Ulate *et al.*, 2015). The latter point is important for this study because we are interested in capturing real-world and not just realistic-looking rainfall patterns at mesoscale resolution. The nudging time-scale was set uniformly to 1 hr wherever nudging was applied.

With the help of a first-order ordinary differential equation, it is straightforward to demonstrate that nudging is effective in constraining the model's seasonally varying climatology and sub-seasonal variations like the MJO. At synoptic time-scales (about 0.5 to 5 days), nudging exerts no coherent impact on tropospheric water vapour mixing ratio because the reference CFSR field and WRF model forcing are not varying in dynamically consistent ways, despite varying much more slowly than the 1 hr time-scale of the nudging process. (Nudging can be loosely compared to chasing after a fluttering butterfly: only the long-time tendency can be followed effectively.) On the time-scale of minutes in moist convective and cloud processes (and thus on the associated spatial scales), the model is practically free to vary in the mid-troposphere despite being nudged. The supplementary material of Fonseca *et al.* (2019) shows that the nudging tendency of water vapour mixing ratio is globally much smaller than the largest term in the precipitable water budget. But this nudging tendency is necessary to correct for the small positive bias in evaporation from subtropical ocean surfaces in the model climatology. The resulting model tropical rainfall is well-verified with observations on pentad time series (Fonseca *et al.*, 2019).

Grid-resolved cloud processes were based on the Goddard six-class microphysics scheme (Tao *et al.*, 1989).

Subgrid-scale cumulus clouds were parametrized using the Betts–Miller–Janjić scheme (Janjić, 1994) but adapted for the tropics (Fonseca *et al.*, 2015), including a precipitating convective cloud scheme to account for the cumulus cloud–radiation feedback (Koh and Fonseca, 2016). The Rapid Radiative Transfer Model for Global (RRTMG) model was used for computing the short-wave and long-wave radiation (Iacono *et al.*, 2008), and the Yonsei University PBL scheme (Hong *et al.*, 2006) was used with Monin–Obukhov surface layer parametrization (Monin and Obukhov, 1954). An interactive prognostic scheme for the sea-surface skin temperature (SSKT) based on Zeng and Beljaars (2005) was used in the model to capture the diurnal variation of the SSKT and allowed its feedback to the air temperature, moisture and radiation fields. The lower boundary condition to the SSKT scheme came from the 6-hourly SST data from CFSR, linearly interpolated in time in order to have a continuously varying thermal constraint on the ocean skin layer.

We are interested in estimating the scaling exponent of the distributions for the cluster area ($A$) and cluster total rain rate ($R$). In a given three-hourly time-slice of the grid-ed model rain product, a single rain cluster was identified as a set of rainy grid pixels with area of $\Delta A$ (the representative area of a model grid point, which was roughly constant about 36 × 36 km) sharing a common border with at least one of its nearest-neighbour rainy grid pixels. The rain area, $A$ in units of $\Delta A$, is the count of rainy grid pixels of the cluster. The total rain rate, $R$ in units of $\Delta A \text{ mm} \cdot \text{hr}^{-1}$ (roughly 360 $\text{m}^3 \cdot \text{s}^{-1}$), is the sum of the rain intensity at all the grid pixels of the cluster. The procedure in data treatment to obtain the distributions of $A$ and $R$ follows closely with that detailed in T17. We refer the readers to Appendix A and Figure A1 for the details of the procedure involved and only highlight here any procedural modifications or departures from T17 when applied to the present analysis.

Four square focus regions of side $Z = 96$ grid pixels over the Indian Ocean (IO) and the western Pacific (WP) warm pools as well as the intertropical convergence zones (ITCZ) in the eastern Pacific (EP) and tropical Atlantic (AO) were identified for analyses (Figure 1). These regions containing deep convection were selected to match as closely as possible the focused regions used in T17 (see also Figure A1). Since the rain clusters found in these regions usually can last for more than 3 hr, many of the rain clusters within a domain for two consecutive three-hour time-slices were not independent. Following T17, we estimated the de-correlation time for the rain clusters for IO, WP, EP and AO to be 9, 6, 8 and 5 days.

Each of the 96 × 96 domains covering the focus regions was further divided into non-overlapping domains of side $Z = 48, 32$ and $24$ respectively (Figure 1b–d) to account for the dependence of the scaling exponents on the domain
<figure><img src="image_6.png" /></figure>

size. As land and coastal regions have the potential of organising rain systems at spatial scales that are determined by the geometry of the coastline and topography, we excluded domains with more than 10% of the area occupied by land.

As our primary interest is in the scaling characteristics of rain clusters arising from convective self-organisation, not all the remaining oceanic domains were used for analysis: some regions like the AO and EP have large swaths of area in which rain arises predominantly from marine stratiform clouds. In our work, we admitted a Z × Z domain for subsequent analyses only if

$$Q(Z; \varepsilon) > Q_0, \qquad (2)$$

where $Q_0$ is a prescribed positive number less than unity (see below) and the statistic Q is a measure of the coverage of convective rain for the domain, defined as

$$Q(Z; \varepsilon) \stackrel{\text{def}}{=} \frac{1}{NZ^2} \sum_{i=1}^{N} \sum_{j=1}^{Z^2} \mathcal{H}(p_{ij} - \varepsilon) \qquad (3)$$

where N is the number of 3-hourly time-slices when there is rainfall somewhere in the domain; $p_{ij}$ is the ratio of rainfall from the convective parametrization to the total rainfall both from the convective and microphysics parametrizations at the jth grid point in the ith time slice; $\varepsilon$ is a prescribed threshold convective fraction of rain (a positive number less than unity); and $\mathcal{H}(k)=1$ if $k>0$ or $\mathcal{H}(k)=0$ otherwise. Q takes into account both the relative contribution of convective (parametrized) rainfall and the frequency and areal coverage by convective rainfall in a domain. The frequency distribution of Q values encountered across all ocean basins and domain sizes for a given $\varepsilon$ (not shown) were found to be bimodal with a minimum around $Q=0.4$. The peak of the distribution at $Q<0.4$ was related to the predominance of shallow convection over EP and AO, while the peak of the distribution at $Q>0.4$ was due to deep convection clusters. The frequency distributions of Q were also fairly robust for a broad range of values of $\varepsilon$ tested. We therefore adopted values of 0.4 and 0.8 for $Q_0$ and $\varepsilon$ respectively.

Figure 1 shows the domains accepted and rejected for the various values of Z. The 96 × 96 domains for EP and AO were excluded because the convective regions are confined within a narrow latitudinal band where the model ITCZ was located (cf. the accepted domains for Z = 48, 32, 24). As the original four 48 × 48 AO domains were rejected by our criteria due to the ITCZ being split almost equally among them, we added two additional 48 × 48 domains containing the ITCZ (heavy solid box in Figure 1b for Z = 48) that satisfy our selection criteria.

We partitioned the cluster area and rain data for all focus regions into 40 non-overlapping temporal data subsets each containing the cluster information for a fixed Coordinated Universal Time (0000 UTC, 0300 UTC, ..., 2100 UTC) with a sampling interval of their respective cluster de-correlation times already mentioned above. The major differences in the details of our analysis method applied here compared with those in T17 have been summarised in Appendix A.

## 3 | PROBABILITY DISTRIBUTIONS

For each focus region, the ensemble mean probability distributions of A and R for all domains of a fixed size Z × Z were obtained separately in each of the 40 non-overlapping temporal data subsets. This ensemble comprised only those Z × Z domains that passed the selection criteria described in the preceding section (Figure 1). Thus, the ensemble members can be considered independent realisations of a true probability distribution, although the 96 × 96 domains for IO and WP form only single-member ensembles.

Figure 2 shows the average of the 40 (correlated) ensemble mean probability distributions for A and R for WP for different Z. Unlike the observed cluster area and rain distributions reported in the literature, which monotonically decrease with increasing A and R (cf. Figure 3 in T17), the model cluster distributions show obvious “bumps” or rapid increase at large A and R. This feature was present for the A and R distributions in the four focus regions across the global tropics. These bumps or sharp increases in the probabilities were the result of the model’s tendency to form a continuous zonal rain band spanning across the domains including the largest 96 × 96 domains which are approximately 30° in extent. This suggests that tropical deep convection in the model tends to be overly organised by planetary-scale processes most of the time. Therefore, for all domain sides Z in the subsequent analysis, we consistently excluded clusters that span across either the east-west sides or north-south sides of a domain, thereby assuring ourselves that we are studying rain clusters that were likely to be self-organised through local processes in the model. The distributions for A and R with these spanning clusters excluded (Figure 3) are monotonic (for all a and for r > 0.1) and the distributions decay faster than power-law behaviour at large values of A and R, in line with observations (e.g. figure 3 of T17). Henceforth, all reference to the cluster distributions refer to the mean ensemble distributions obtained by excluding any east-west domain-spanning clusters.
<figure><img src="image_7.png" /></figure>

<figure><img src="image_8.png" /></figure>

**FIGURE 2** Double logarithmic plots of the mean of the ensemble probability distribution of (a) rain cluster area *A*, and (b) total cluster rain rate *R*, over WP of Figure 1. *A* and *R* are measured in units of Δ*A* and Δ*A* mm·hr⁻¹ respectively, where Δ*A* is the (assumed constant) representative area of an evaluation grid point in the model output. The domain size *Z* × *Z* is specified in the legend

<figure><img src="image_9.png" /></figure>

<figure><img src="image_10.png" /></figure>

**FIGURE 3** (a,b) Same as Figure 2 but with the east-west and north-south spanning clusters excluded when computing the statistics (see the main text for details). The insert at (b) is the magnified view that emphasises the dependence of the scaling exponent on the domain size. The grey dotted lines in all the plots are reference power law with their respective scale exponents, ζ<sub>S</sub>, estimated for the *Z* = 96 domain for WP (see the main text). These reference lines are drawn as a guidance as to where the simple scaling regime are approximately assumed in this study

From Figure 3a, the distributions for *A* for all domain sizes appear to follow simple power-law scaling beyond *a* ~ 10Δ*A*, though we note that the scaling regime even for the largest domain size (Z = 96) seemed to span at most only one order of magnitude. The distribution for *R* is more complicated. For very small values of *R* (≤10⁻⁴ Δ*A* mm·hr⁻¹), the distribution is highly variable (Figure 3b). It was found that most of these light rain events were associated with clusters at the limit of model resolution (*a* = 1), and the high variability in the *R* distributions is due to the discretization of the very light rain produced by the cloud microphysics and the cumulus parametrizations.

Unlike the case for *A*, it is more difficult to discern a power-law regime visually for the distribution for larger *R*; the steepening in the gradient of the distribution around *r*~0.5 Δ*A* mm·hr⁻¹ is gradual and a power-law regime is not easily discerned for domains with *Z* < 96 from the plot (Figure 3b). However, we proceed to assume a simple scaling ansatz for the distribution of *R* (for large enough *R*)

since it is possible that the magnitude of the upper cut-off, $r_c$, where the finite system size starts deviating from the strict power law is relatively small (compared to $r_0$, the lower cut-off where simple scaling ansätze are valid; see Equation 1) making it difficult to discern a clearly defined scaling regime within the available domain size. Moreover, the observed distribution of *R* over the same oceanic domains exhibited simple scaling (T17) and hence it is not unreasonable to use simple scaling for *R a priori*.

# 4 | SCALING EXPONENTS

## 4.1 | Exponent for the probability distributions

For each focus region, the scaling exponents for both the *A* and *R* distributions for the different *Z* × *Z* domains were estimated by fitting the cluster ensemble distribution to a
<figure><img src="image_11.png" /></figure>

<figure><img src="image_12.png" /><figcaption>FIGURE 4 (a–d) Robust iterative least-square linear regression with Huber weights of scaling exponents versus 1/Z for the different focus regions for rain cluster size A (black open squares), and total cluster rain rate R (grey crosses). Each data point in the plots is the ensemble-averaged scaling exponent, $\tau_S(Z, S_0)$, S representing A or R. (e) Estimate for $\zeta_A$ (black) and $\zeta_R$ (grey) for individual focus regions. These point estimates are the regression intercepts of (a–d). The error bars represent the limits at 95% confidence level of the estimates</figcaption></figure>

distribution model:

$$f(s) = C_S s^{-\tau_S} \exp \left[ \left( \frac{s}{S_c} \right)^{\nu_S} \right], s \ge s_0; S \in \{A, R\}, \quad (4)$$

where $f(s)$ is the (discrete) probability for the case when S is A (measured in $\Delta A$), or the probability density function for the case when S is R (measured in $\Delta A \text{ mm} \cdot \text{hr}^{-1}$); $C_S$ is a normalisation constant; and $\tau_S$, $\nu_S$ and $s_c$ are distribution parameters determined using the Method of Maximum Likelihood for a prescribed $s_0$, labelled as $a_0$ and $r_0$ respectively for A and R. In this work, we selected six $a_0$ approximately separated by regular logarithmic intervals from 10 to 50. Similarly, six $r_0$ at regular logarithmic separations from 0.5 to $5 \Delta A \text{ mm} \cdot \text{hr}^{-1}$ were selected. We obtained one estimate of the scaling exponent, $\tau_S(Z; s_0)$, $S \in \{A, R\}$, from each of the 40 ensemble distributions computed from the cluster data subsets (see above). The best point estimate of $\tau_S$ for each Z and $s_0$ was then obtained from these 40 $\tau_S(Z; s_0)$ estimates following T17.

Visual inspection of the plots of the scaling exponents versus the reciprocal of the domain size Z, suggested that there is a dependence of the scaling exponents on the domain size in which the clusters were identified, for IO, WP, EP (Figure 4a–c). Such a dependence can also be

discerned from the slopes of the near-linear segments in the plots of the probability distributions of A and R for the three focused regions (cf. Figure 3). Following T17, we estimated $\zeta_S$ as the thermodynamic limit of $\tau_S$ (i.e. $\zeta_S = \lim_{Z \to \infty} \tau_S$) using the intercept of a robust linear regression with Huber weights of the point estimates of $\tau_S(Z; S_0)$ with 1/Z as in T17 (Figure 4a–d). All the regressions for the focus regions were significant (p-value < .05) except for AO. Figure 4e shows the estimated $\zeta_S$, and their values are tabulated in Table 1. These scaling exponents identified in this study are compared with the corresponding values in the four ocean basins in T17 (Appendix B and Figure B1).

## 4.2 | Scaling exponent of conditional mean of R

For each focus region and each domain size Z, the sample estimates of $E(R|a)$, the conditional mean of the cluster rain amount given the cluster area $A = a$, were calculated from each of the 40 data subsets. For a given Z, the conditional mean was computed from the clusters within those domains that passed the selection criteria described in Section 2. Figure 5a–d show the double logarithmic plots of the $E(R|a)$ versus a for the statistics from one of
<figure><img src="image_13.png" /></figure>

<table><caption>TABLE 1 The estimated scaling exponents</caption><thead><tr><th></th><th>IO</th><th>WP</th><th>EP</th><th>AO</th></tr></thead><tbody><tr><th>ζ<sub>A</sub></th><td>2.12 ± 0.10</td><td>2.48 ± 0.09</td><td>2.37 ± 0.16</td><td>2.12 ± 0.08</td></tr><tr><th>ζ<sub>R</sub></th><td>1.89 ± 0.07</td><td>2.21 ± 0.12</td><td>2.38 ± 0.17</td><td>1.97 ± 0.21</td></tr><tr><th>β</th><td>1.31 ± 0.08</td><td>1.28 ± 0.04</td><td>1.33 ± 0.07</td><td>1.42 ± 0.54</td></tr><tr><th>β̂</th><td>1.26 ± 0.18</td><td>1.22 ± 0.18</td><td>1.00 ± 0.24</td><td>1.16 ± 0.34</td></tr></tbody></table>

Note: Bold values in the first two rows are estimates for which the regression p-values are less than .05. β is estimate of exponent b in Equation (5). The last row shows the estimates of β using the scaling relation Equation (6) in the main text.

<figure><img src="image_14.png" /><figcaption>FIGURE 5 (a–d) The conditional average of cluster total rain rate with different cluster area for the various domain sizes over IO. The solid line shows the linear regression. (e) Estimate β for the scaling exponents of individual focus regions. The horizontal lines show the range of values where the separate estimates for all focus regions overlap. The error bars represent the limits at 95% confidence level of the estimates</figcaption></figure>

the subsets for IO for the different domain sizes, which strongly suggest that $E(R|a)$ has a power-law scaling. Such scaling has also been reported in T17 for the observed oceanic rain clusters. Figure 5 is fairly representative of the $E(R|a)$ for clusters from the remaining data subsets as well as for other focus regions. T17 reported that the observed $E(R|a)$ have two scaling regimes, with rain clusters smaller than length-scale of ~320 km having an exponent of 4/3 and larger clusters having $E(R|a)$ proportional to the rain cluster size. In contrast the $E(R|a)$ from the model did not have any discernible transition in its scaling exponent (cf. Figure 5a–d). Note that a rain cluster of length-scale ~320 km is about 79 $\Delta A$ in our model configuration lying mid-range in Figure 5a–d. Hence, we adopted a simpler

model than T17:

$$E(R|a) = ka^b, \qquad (5)$$

where k and b are positive numbers to be determined. For each focused region, an estimate of the scaling exponent b was obtained from each data subset via a linear regression of $log E(R|a)$ against $log a$. The procedure of estimating b was repeated for all Z and a representative $b(Z)$ was subsequently determined as the sample mean of the 40 separate b estimates. Separate b values were estimated for each Z in anticipation of a functional dependence of b with Z as in the case of $\zeta_s$. There were no $b(Z=96)$ estimates for EP and AO however, as the 96 × 96 domains for both focus regions did not pass the selection criteria described in Section 2.
<figure><img src="image_15.png" /></figure>

However, the p-values of the linear regression of the b(Z) against 1/Z (not shown) as in the preceding section were greater than .05 for all focus regions, providing no significant evidence that the exponent b depends on Z. As such, we estimated the scaling exponent as the sample mean β of the available b(Z) (Figure 5e and Table 1).

As a consistency check, from the scaling relation in T17,

$$ \hat{\beta} = \frac{1 - \zeta_A}{1 - \zeta_R}, \qquad (6) $$

another estimate $\hat{\beta}$ of the scaling exponent for E(R|a) against a was computed for each ocean basin using the currently estimated $\zeta_A$ and $\zeta_R$. Except for the near-miss in EP, $\hat{\beta}$ of all ocean basins practically agree within error bounds with β directly obtained from the scaling of E(R|a) described earlier (see Table 1). This is true because Equation (6) is derived from statistical arguments alone.

Our results for β in Figure 5e are consistent with the view that the scaling exponent is the same within the error bounds across all the ocean basins, although we note the relatively large uncertainty in the estimate for AO. One estimate for a representative value for all the oceans is taking the mean of the individual estimates which gives $\bar{\beta} = 1.36 \pm 0.05$. We include AO in the computation because, notwithstanding the lower than 95% confidence level for $\zeta_A$, $\zeta_R$ and hence $\hat{\beta}$ estimates in AO, this focus region shows the same insensitivity of the b estimates to domain size Z as the other focus regions. Hence, despite the differences in $\zeta_A$ and $\zeta_R$ across the ocean basins, we think that β is the same across all the ocean basins. Furthermore, the estimate of $\bar{\beta}$ is the same within estimation errors as the observed β value of $1.33 \pm 0.03$ reported in T17.

## 5 | FURTHER DISCUSSIONS

One of the requirements for the rain clusters to be universal is that the scaling exponents, $\zeta_A$ and $\zeta_R$, should be identical across all the oceanic basins. Our results show that the modelled rain clusters do not exhibit universality as observed in T17. First, the functional dependence $\zeta_A$ and $\zeta_R$ with the domain size Z is different across the ocean basins: the scaling exponents with domain size for AO appears to be independent of Z in contrast with the other three ocean basins. This result suggests that the scaling behaviour of model rain cluster distribution over the tropical Atlantic is different from those in the Indian and Pacific Oceans. Moreover, there is no one constant value for $\zeta_A$ or $\zeta_R$ that lies within the error bounds for all three ocean basins, IO, WP and EP in Figure 4e. The scaling exponents in each ocean basin from the model are larger than the reported values of $1.66 \pm 0.06$ (~5/3) and

$1.48 \pm 0.13$ (~3/2) for $\zeta_A$ and $\zeta_R$ respectively in T17, and for $\zeta_R$ reported in Quinn and Neelin (2017b) based on Special Sensor Microwave Imager data. The $\zeta_R$ values are also higher than the globally average values (1.39/1.36 simulated by model of resolution 25/50 km) reported in Quinn and Neelin (2017a). Hence, the model reproduced too few larger or heavier rain clusters and too many smaller or lighter rain clusters over the tropical oceans than observed. However, we note that the set of scaling exponents for IO curiously agrees with those reported by Peters et al. (2010) for the western Pacific, perhaps by coincidence.

The E(R|a) for observed rain clusters of length-scale of ~320 km and beyond were reported to scale linearly with the cluster area in T17. This was interpreted by the authors as the tendency for the larger rain clusters to grow without intensification since the average intensity within the cluster, E(R|a)/a, is approximately constant. On the other hand, the clusters smaller than ~320 km tend to grow and intensify simultaneously as β ~ 4/3 in the observations. Although β in the model generally agrees with the observed value in T17 for the smaller rain clusters, no crossover of the scaling exponent for E(R|a) was seen for the larger rain clusters in the model (Figure 5). It therefore appears that the modelled rain clusters continue to intensify as the rain clusters grow in size, but less and less effectively as $0 < (\beta - 1) < 1$. This is so even for superclusters of ~1,000 km in size (~700 ΔA), suggesting that the self-organisation of the large-scale precipitation in the model departs significantly from reality.

Our results also suggest that the model could over-predict the rain intensity for synoptic rain clusters. A possible consequence is that the anomalous diabatic heating by these synoptic rain clusters would reduce the static stability of the atmosphere and consequently reduce the equivalent depth of the convectively coupled equatorial waves in the model. In addition, such large-scale release of extra latent heat would also necessitate an anomalous large-scale ascent within the rain clusters and a compensating large-scale descent outside the rain clusters. We suspect that these anomalous large-scale circulations would potentially confine further convective development to within the clusters, causing them to grow into unrealistic basin-wide continuous rain bands highlighted in Section 2. As we know that energy interacts across scales in nonlinear flow dynamics, the emergence of excessive spectral power at large scales would affect the spectral power at small scales. Consequently, it is hard to expect the scaling exponents $\zeta_A$ and $\zeta_R$ to be well-captured in the model despite having excluded those erroneous basin-wide rain bands from our statistical analyses.

Although there is no common value of $\zeta_A$ or $\zeta_R$ amongst the ocean basins in the model, there is a common β, and its value agrees with observations for mesoscale
rain clusters. In addition, our analysis suggests that perhaps $\beta$ plays a more fundamental role in the sense that it directly concerns the intensification (or dissipation) of rain clusters as they grow (or shrink). We suspect the agreement of $\beta$ between the observations and the model could be because the scaling of $E(R|a)$ with $a$ depends on a small set of dynamical processes of mesoscale deep convection that the model captures adequately as we qualitatively illustrate next.

For a rain cluster of area $a$, whose length-scale, $\ell \equiv a^{1/2}$, is much larger than the scale height of the rain cluster, $H$, it is not unreasonable to view the rain cluster as a prototypical turbulent eddy in a two-dimensional turbulent fluid, where energy is transferred upscale (Kraichnan, 1967). In this turbulent framework, energy is injected within the rain cluster at the convective (small) scale as a result of the deep convection embedded within the rain clusters. The energy generated is transferred upscale to the larger-scale flow within the rain cluster through the nonlinearity of the turbulent flow. This view of inverse energy cascade from convective scale to the larger horizontal flow structure in the rain cluster is consistent with the results from past investigations (Schmidt and Cotton, 1990; Moncrieff, 1992; Pandya and Durran, 1996) which show that the circulation within the mesoscale convective system (MCS; see Figure 6b) is a gravity wave response of a latent heating profile of the MCS with a top-heavy heating profile in the vertical. For a rain cluster with cluster rain amount $R$ (unit of volume per time), we assert that this rate of energy transferred per unit mass of air, $\epsilon$, across the spatial scales is equal to the rate of latent heat released per unit mass of air within the rain cluster:

$$ \epsilon = \frac{\rho_{\ell} L_{\nu} r}{\bar{\rho} A H}, \qquad (7) $$

where $\rho_{\ell}$ and $\bar{\rho}$ are the density of liquid water and a representative air density within the cluster respectively, and $L_{\nu}$ is the latent heat of vaporisation of water. If the rain clusters are turbulent eddies within the inertial subrange of a homogeneous turbulent fluid, the kinetic energy per unit mass is characterised only by $\ell (= a^{1/2})$ and $\epsilon$:

$$ u^2 \sim \epsilon^{2/3} a^{1/3}, \qquad (8) $$

where $u$ is the characteristic horizontal flow speed within the rain cluster. Using $r \propto a^{\beta}$ and substituting $\epsilon$ in Equation (8) with Equation (7), and ignoring all the constants (which we assume all variables other than $u$, $r$ and $a$ are), we have a relationship between the (horizontal) kinetic energy and the cluster area:

$$ u^2 \sim a^{\frac{2}{3}(\beta-1)+\frac{1}{3}}. \qquad (9) $$

<figure><img src="image_16.png" /><figcaption>FIGURE 6 (a) Idealisation of a top view of an MCS divided into convective (dark grey) and stratiform (hatched) regions. In the discussion in the main text, it is assumed that the convective region has a half-perimeter, $p_c$, and is embedded with $n$ convective centres with strong updraughts of length-scale $D$. (b) Schematic diagram showing the airflow relative to a two-dimensional, steady-state mesoscale convective system in a large-scale environment of given wind shear. The environmental air entering the updraught is potentially unstable, and there is a pressure and temperature decrease across the system at the upper level going from A to B. The streamlines are those required by conservation of mass, momentum, entropy and vorticity. Adopted from Houze (2004)</figcaption></figure>

Alternatively, we may appeal to the archetypal MCS shown in Figure 6 as a model for the rain cluster to obtain a scaling relation between kinetic energy and cluster area. Referring to the schematics of Figure 6b, the upper-level outflow of the MCS is maintained by the pressure difference, $\Delta p$, between the region within the MCS (A in Figure 6b) and the environmental air outside the rain cluster (B in Figure 6b). Using the Bernoulli equation:

$$ \frac{1}{2} u^2 \approx \frac{1}{\rho} \Delta p. \qquad (10) $$

Assuming that this pressure difference is primarily due to the diabatic heating of the deep convection within the rain cluster:

$$ \frac{1}{\rho} \Delta p \approx c_p \Delta T \approx \epsilon \tau, \qquad (11) $$

where $c_p$ is the specific heat capacity per unit mass of the air at constant pressure; $\Delta T$ the horizontal temperature difference between A and B in Figure 6b; $\tau$, a characteristic time-scale of the MCS. Adopting an argument put forth by Emmanuel (1994), we assumed that in time $\tau$, the deep convection within the MCS would exhaust conditionally unstable air near the surface of thickness $d$, covering an area similar to (or at least proportional to) the rain cluster area:

$$ \tau \approx \frac{ad}{p_c v_c} \qquad (12) $$

where $v_c$ is the characteristic horizontal air speed for the thickness $d$ feeding into the convective area and $p_c$ is the perimeter of the convective “patch” (see Figure 6a). In the absence of any background wind, that is, in a purely
<figure><img src="image_17.png" /></figure>

convective overturning case, $v_c$ would scale as the vertical updraught, which in turn scales approximately as the square root of the convectively available potential energy (CAPE). We make two further assumptions: (a) the CAPE does not vary much over tropical oceans (Emanuel, 1994), and (b) $p_c \sim \ell^K (= a^{K/2})$, for some dimension exponent K. Several observational studies (Cahalan and Joseph, 1989; Mapes and Houze, 1993; Benner and Curry, 1998; Peters and Neelin, 2006) reported fractal dimension for the perimeter of clouds and rain clusters, with K ranging from 1.37 to 1.5. Taking all the assumptions described above together, we obtained from Equations (10)–(12):

$$u^2 \sim a^{\beta-K/2}. \qquad (13)$$

Equating the exponents of a in Equations (9) and (13), we have

$$\beta = \frac{3}{2}K - 1. \qquad (14)$$

Hence $\beta \ge 1$ if $K \ge 4/3$. For $\beta = 4/3$ as revealed by the TRMM observations analysed in T17 and captured by the WRF model here, $K = 14/9$ (~1.56), which is just above the range of the reported values. It is interesting also to note Peters and Neelin (2006) argued that the deep convection clusters over the tropical oceans could be a gradient percolation, and consequently $K = 4/3$, which gives $\beta = 1$ according to Equation (14).

In short, the phenomenological argument above explains why $\beta$ can be greater than unity, and the value of $\beta$ depends on the (fractal) dimension of the perimeter of the organised deep convection within the rain cluster.

## 6 | CONCLUSION

The scaling exponents of the area and total rain amount for tropical oceanic rain clusters in the predominantly deep convective regions from a numerical weather model were obtained. The model oceanic clusters did not exhibit the universality reported in T17 for observed rain clusters: the scaling exponents $\zeta_A$ and $\zeta_R$ do not agree across the ocean basins unlike what is observed. In addition, the estimated values of both $\zeta_A$ and $\zeta_R$ from all the ocean basins were larger than their corresponding observed values reported in T17. The model also produced unrealistic basin-wide clusters within the ITCZs in all the ocean basins, which were excluded from our estimates of scaling exponents.

The scaling of the model conditional mean rainfall, $E(R|a)$, with the rain cluster area $a$ was better replicated. A robust constant value was found for the scaling exponent $\beta$ across all the ocean basins in the model

and its value $1.36 \pm 0.05$ is close to the observed $\beta$ of about 4/3 for mesoscale clusters of sizes about 36 km to 320 km. Hence, the model was able to reproduce the correct increase in the average rain intensity with increasing rain cluster area. We interpret this as an indication that the model has captured the fundamental and universal nature of mesoscale convective flow dynamics, as illustrated conceptually by a dimensional argument based on a simple model.

Notwithstanding the above success, the model failed to reproduce the crossover in the scaling of $E(R|a)$ as cluster size increases beyond ~320 km as seen in T17, where the average rain intensity of the synoptic-scale rain clusters is observed to be independent of the cluster area (i.e. $\beta \sim 1$). In the model, the intensification of rain with cluster size continues with the same $\beta \sim 4/3$ for the synoptic rain clusters. We think that this may be the reason why unrealistic basin-wide rain bands form in the model, leading indirectly to the model's scaling exponents $\zeta_A$ and $\zeta_R$ deviating from reality.

The modelling community has advanced rapidly in recent years, partly in association with the activities in the sixth generation of the Climate Model Intercomparison Project (CMIP6). For example, many models in the High-Resolution Model Intercomparison Project (High-ResMIP: Haarsma et al., 2016) have horizontal resolution of the order 20–50 km, which is comparable to that in the WRF applied in this study. Most of these high-resolution models still apply convection parametrization, thus the simulated rainfall is the sum of that explicitly resolved and that from the parametrization scheme is similar, as in our WRF model. Our theoretical arguments have suggested that the scaling exponent $\beta = 4/3$ (Equation 5) depended only on a few fundamental features of the known mesoscale convection system, hence we expect most models to reproduce this exponent. It would be important to see whether these models also produce super (domain-spanning) clusters as what we observed with the present model, and whether they are able to reproduce the observed scaling transition to the $\beta = 1$ scaling (i.e. simple merging process). If there are a variety of different scaling behaviours across these models, further investigations could be followed up to understand how much of the differences in the modelled scaling properties of the rain clusters may be attributed to the model details, for example the convection parametrization used.

Advances in climate modelling include the global cloud-resolving models such as those applied in the Dynamics of the Atmospheric general circulation Modelled On Non-hydrostatic Domains (DYAMOND: Stevens et al., 2019; see also Satoh et al., 2019). As all the cloud microphysical processes are modelled explicitly, the scaling behaviours of the rain clusters from these
models may be used to assess the impact of cumulus parametrization on the scaling behaviour of the tropical rain clusters. Our suggested scaling exponents for the tropical rain clusters may also be useful in assessing idealised Radiative-Convective Equilibrium (RCE) models, such as those model intercomparison projects in the RCE-MIP (Wing et al., 2018). In the RCE-MIP, cloud cluster aggregation was identified as a feature to be assessed using metrics including organisation index, subsidence fraction and autocorrelation length (Wing et al., 2018). We believe the rain cluster scaling exponents of A and R can be considered as another set of metrics to evaluate cloud clustering and aggregation in these RCE models.

Although extreme precipitation is not the focus in this study, our findings have implications for the recent studies on the relationship between extreme precipitation and convective organisation/aggregation (e.g. Semie and Bony, 2020). One of the debates in these studies is whether precipitation extreme increases with convective organisation with the definition of extreme varying from instantaneous to short time-scale (e.g. daily), and spatially considering an area-specific extreme or only that within deep convection (see review in Pendergrass, 2020). So far the result of $\beta \sim 4/3$ from both observations and our current model analysis suggests that the extreme part of the within-cluster rainfall intensity ($R/A \sim A^{1/3}$) likely would increase with the degree of organisation, at least for those mesoscale convective clusters or smaller. However, a more explicit comparison between our identified organised clusters with other measures of organisation is needed. On the other hand, the scale regime of $E(R|a)$ changes to $\beta = 1$, at least for the larger observed rain clusters beyond length-scale of $\sim 300$ km (T17). If the cluster rain within this regime grows in size mainly through agglomeration of smaller rain clusters, the domain-averaged extreme precipitation intensity ($R/A \sim \text{const}$) might very well be independent of the fractional area of rain. However we cautioned that the above speculations of what our analysis of the scaling exponent of first-order moment of the (conditional) rain cluster distribution ($E(R|a)$) would reveal about the relation between extreme rainfall and the organisation of a rain cluster had implicitly assumed that the higher-order moments of the rain cluster rainfall distributions were independent, or at the very least, not a strong function of the cluster area. Further investigations would be needed to determine whether these assumptions are valid.

## ACKNOWLEDGEMENTS

This research is in collaboration with the Center for Environmental Sensing and Modeling (CENSAM) and is supported by the National Research Foundation Singapore under its Campus for Research Excellence and

Technological Enterprise programme. CENSAM is an interdisciplinary research group of the Singapore-MIT Alliance for Research and Technology. We would like to thank the two anonymous reviewers whose comments have improved the quality of our manuscript substantially.

## ORCID

Kevin K. W. Cheung https://orcid.org/0000-0002-1204-3089

## REFERENCES

Ahmed, F. and Neelin, J.D. (2019) Explaining scales and statistics of tropical precipitation clusters with a stochastic model. *Journal of Climate*, 76, 3063–3087.

Andrade, R.F.S., Schellnhuber, H.J. and Claussen, M. (1998) Analysis of rainfall records: possible relation to self-organized criticality. *Physica A*, 354, 557–568.

Andrews, D.G., Holton, J.R. and Leovy, C.B. (1987) *Middle Atmosphere Dynamics*. New York, NY: Academic Press.

Bak, P., Tang, C. and Wiesenfeld, K. (1987) Self-organized criticality: an explanation of the 1/f noise. *Physical Review Letters*, 59(4), 381–384.

Benner, T.C. and Curry, J.A. (1998) Characteristics of small tropical cumulus clouds and their impact on the environment. *Journal of Geophysical Research*, 103(D22), 28753–28767. https://doi.org/10.1029/98JD02579.

Biello, J.A. and Majda, A.J. (2005) A new multiscale model for the Madden–Julian Oscillation. *Journal of the Atmospheric Sciences*, 62, 1694–1720.

Cahalan, R.F. and Joseph, J.H. (1989) Fractal statistics of cloud fields. *Monthly Weather Review*, 117, 261–272.

Christensen, K. and Moloney, N.R. (2005) *Complexity and Criticality*. London: Imperial College Press.

Emanuel, K.A. (1994) *Atmospheric Convection*. New York, NY: Oxford University Press.

Fang, X. and Kuo, Y.H. (2015) A new generic method for quantifying the scale predictability of the fractal atmosphere: applications to model verification. *Journal of the Atmospheric Sciences*, 72, 1667–1688.

Fonseca, R., Koh, T.Y. and Teo, C.K. (2019) Multi-scale interactions in a high-resolution tropical-belt experiment and observations. *Climate Dynamics*, 52(5), 3503–3532. https://doi.org/10.1007/s00382-018-4332-y.

Fonseca, R.M., Zhang, T. and Koh, T.Y. (2015) Improved simulation of precipitation in the tropics in WRF using a modified BMJ scheme. *Geoscientific Model Development*, 8, 2915–2928.

Haarsma, R., Roberts, M.J., Luigi Vidale, P., Senior, C.A., Bellucci, A., Bao, Q., Chang, P., Corti, S., Fuckar, N.S., Guemas, V., von Hardenberg, J., Hazeleger, W., Kodama, C., Koenigk, T., Leung, R., Lu, J., Luo, J.J., Mao, J., Mizielinski, M.S., Mizuta, R., Nobre, P., Satoh, M., Scoccimarro, E., Semmier, T., Small, J. and von Storch, J.S. (2016) High-resolution model intercomparison project (High-ResMIP v1.0) for CMIP6. *Geoscientific Model Development*, 9, 4185–4208.

Hong, S.Y., Noh, Y. and Dudhia, J. (2006) A new vertical diffusion package with an explicit treatment of entrainment processes. *Monthly Weather Review*, 134, 2318–2341.
Hottovy, S.A. and Stechmann, S.N. (2015) A spatiotemporal stochastic model for tropical precipitation and water vapor dynamics. *Journal of the Atmospheric Sciences*, 72, 4721–4738.

Houze, R.A., Jr. (2004) Mesoscale convective systems. *Reviews of Geophysics*, 42, RG4003. https://doi.org/10.1029/2004RG000150.

Iacono, M.J., Delamere, J.S., Mlawer, E.J., Shephard, M.W., Clough, S.A. and Collins, W.D. (2008) Radiative forcing by long-lived greenhouse gases: calculations with the AER radiative transfer models. *Journal of Geophysical Research*, 113(D13), D13103. https://doi.org/10.1029/2008JD009944.

Janjić, Z.I. (1994) The step-mountain eta coordinate model: further developments of the convection, viscous sublayer and turbulence closure schemes. *Monthly Weather Review*, 122, 927–945.

Kahn, B.H., Teixeira, J., Fetzer, E.J., Gettelman, A., Hristova-Veleva, S.M., Huang, X., Kochanski, A.K., Mohler, M., Krueger, S.K., Wood, R. and Zhao, M. (2011) Temperature and water vapor variance scaling in global models: comparisons to satellite and aircraft data. *Journal of the Atmospheric Sciences*, 68, 2156–2168.

Kiladis, N., Wheeler, M.C., Haertel, P.T., Straub, K.H. and Roundy, P.E. (2009) Convectively coupled equatorial waves. *Reviews of Geophysics*, 47(2), RG2003. https://doi.org/10.1029/2008RG000266.

Koh, T.-Y. and Fonseca, R. (2016) Subgrid-scale cloud-radiation feedback for the Betts-Miller-Janjić convection scheme. *Quarterly Journal of the Royal Meteorological Society*, 142(695), 989–1006.

Kraichnan, R. (1967) Inertial ranges in two-dimensional turbulence. *Physics of Fluids*, 10, 1417–1423.

Mapes, B.E. and Houze Jr, R.A. (1993) Cloud clusters and superclusters over the oceanic warm pool. *Monthly Weather Review*, 121, 1398–1415.

Moncrieff, M.W. (1992) Organized convective systems: archetypal dynamical models, mass and momentum flux theory, and parametrization. *Quarterly Journal of the Royal Meteorological Society*, 118(507), 819–850.

Monin, A.S. and Obukhov, A.M. (1954) *Basic laws of turbulent mixing in the ground layer of the atmosphere*. Akademiia nauk SSSR, Vol. 151. Boston, MA: Geofizicheskii institut, pp. 163–187.

Nakazawa, T. (1988) Tropical super clusters within intraseasonal variations over the western Pacific. *Journal of the Meteorological Society of Japan*, 66, 823–839.

Neelin, J.D., Peters, O., Lin, J.W.B., Hales, K. and Holloway, C.E. (2008) Rethinking convective quasi-equilibrium observational constraints for stochastic convective schemes in climate models. *Philosophical Transactions of the Royal Society A*, 366, 2581–2604.

O'Brien, T.A., Li, R., Collins, W.D., Rauscher, S.A., Ringler, T.D., Taylor, M., Hagos, S.M. and Leung, R. (2013) Observed scaling in clouds and precipitation and scale incognizance in regional to global atmospheric models. *Journal of Climate*, 26, 9313–9333.

Pandya, R. and Durran, D. (1996) The influence of convectively generated thermal forcing on the mesoscale circulation around squall lines. *Journal of the Atmospheric Sciences*, 53, 2924–2951.

Pendergrass, A.G. (2020) Changing degree of convective organization as a mechanism for dynamic changes in extreme precipitation. *Current Climate Change Reports*, 6, 47–54.

Peters, O., Deluca, A., Corral, Á., Neelin, J.D. and Holloway, C.E. (2010) Universality of rain event size distributions. *Journal of Statistical Mechanics: Theory and Experiment*, 11, P11030.

Peters, O. and Neelin, J.D. (2006) Critical phenomena in atmospheric precipitation. *Nature Physics*, 2(6), 393–396.

Peters, O., Neelin, J.D. and Nesbitt, S.W. (2009) Mesoscale convective systems and critical clusters. *Journal of the Atmospheric Sciences*, 66, 2913–2924.

Quinn, K. and Neelin, J.D. (2017a) Distribution of tropical precipitation cluster power and their changes under global warming. Part I: Observational baseline and comparison to a high-resolution atmospheric model. *Journal of Climate*, 30, 8033–8044.

Quinn, K. and Neelin, J.D. (2017b) Distribution of tropical precipitation cluster power and their changes under global warming. Part II: Long-term time dependence in coupled model intercomparison project phase 5 models. *Journal of Climate*, 30, 8045–8059.

Saha, S., Moorthi, S., Pan, H.-L., Wu, X.R., Wang, J.D., Nadiga, S., Tripp, P., Kistler, R., Woollen, J., Behringer, D., Liu, H.X., Stokes, D., Grumbine, R., Gayno, G., Wang, J., Hou, Y.-T., Chuang, H.-Y., Juang, H.-M.H., Sela, J., Iredell, M., Treadon, R., Kleist, D., Van Delst, P., Keyser, D., Derber, J., Ek, M., Meng, J., Wei, H.L., Yang, R.Q., Lord, S., van den Dool, H., Kumar, A., Wang, W.Q., Long, C., Chelliah, M., Xue, Y., Huang, B.Y., Schemm, J.-K., Ebisuzaki, W., Lin, R., Xie, P.P., Chen, M.Y., Zhou, S.T., Higgins, W., Zou, C.-Z., Liu, Q.H., Chen, Y., Han, Y., Cucurull, L., Reynolds, R.W., Rutledge, G. and Goldberg, M. (2010) The NCEP climate forecast system reanalysis. *Bulletin of the American Meteorological Society*, 91, 1015–1057.

Satoh, M., Stevens, B., Judt, F., Khairoutdinov, M., Lin, S.J., Putman, W.M. and Duben, P. (2019) Global cloud-resolving models. *Current Climate Change Reports*, 5, 172–184.

Schmidt, J.M. and Cotton, W.R. (1990) Interactions between upper and lower tropospheric gravity waves on squall line structure and maintenance. *Journal of the Atmospheric Sciences*, 47, 1205–1222.

Semie, A.G. and Bony, S. (2020) Relationship between precipitation extremes and convective organization inferred from satellite observations. *Geophysical Research Letters*, 47, e2019GL086927. https://doi.org/10.1029/2019GL086927.

Skamarock, W.C. (2004) Evaluating mesoscale NWP models using kinetic energy spectra. *Monthly Weather Review*, 132, 3019–3032.

Skamarock, W.C., Klemp, J.B., Dudhia, J., Gill, D.O., Barker, D.M., Duda, M.G., Huang, X.-Y., Wang, W. and Powers, J.G. (2008) *A description of the Advanced Research WRF Version 3*. NCAR Technical Note TN-4175+STR.

Stechmann, S.N. and Neelin, J.D. (2014) First-passage-time prototypes for precipitation statistics. *Journal of the Atmospheric Sciences*, 71, 3269–3291.

Stevens, B., Satoh, M., Auger, L., Biercamp, J., Bretherton, C.S., Chen, X., Duben, P., Just, F., Khairoutdinov, M., Klocke, D., Dodama, C., Kornblueh, L., Lin, S.J., Neumann, P., Putman, W.M., Rober, N., Shibuya, R., Vanniere, B., Vidale, P.L., Wedi, N. and Zhou, L. (2019) DYAMOND: the dynamics of the atmospheric general circulation modeled on non-hydrostatic domains. *Progress in Earth and Planetary Science*, 6, 61. https://doi.org/10.1186/s40645-019-0304-z.

Tao, W.K., Simpson, J. and McCumber, M. (1989) An ice-water saturation adjustment. *Monthly Weather Review*, 117, 231–235.

Teo, C.K., Huynh, H.N., Koh, T.Y., Cheung, K.K.W., Legras, B., Chew, L.Y. and Norford, L. (2017) The universal scaling characteristics of tropical oceanic rain clusters. *Journal of Geophysical Research: Atmospheres*, 122, 5582–5599. https://doi.org/10.1002/2016JD025921.

Traxl, D., Boers, N., Rheinwalt, A., Goswami, B. and Kurths, J. (2016) The size distribution of spatiotemporal extreme clusters around
the globe. *Geophysical Research Letters*, 43, 9939-9947. https://doi.org/10.1002/2016GL070692.

Tulich, S.N. and Mapes, B.E. (2008) Multiscale convective wave disturbances in the tropics: insights from a two-dimensional cloud-resolving model. *Journal of the Atmospheric Sciences*, 65, 140–155.

Ulate, M., Zhang, C. and Dudhia, J. (2015) Role of water vapour and convection-circulation decoupling in MJO simulations by a tropical channel model. *Journal of Advances in Modeling Earth Systems*, 7, 692–711.

Wing, A.A., Reed, K.A., Satoh, M., Stevens, B., Bony, S. and Ohno, T. (2018) Radiative-convective equilibrium model intercomparison project. *Geoscientific Model Development*, 11, 793–813.

Zeng, X. and Beljaars, A.C.M. (2005) A prognostic scheme of sea surface skin temperature for modeling and data assimilation. *Geophysical Research Letters*, 32(14), L14605. https://doi.org/10.1029/2005GL023030.

**How to cite this article:** Teo C-K, Koh T-Y, Cheung KKW, *et al*. Scaling characteristics of modelled tropical oceanic rain clusters. *QJR Meteorol Soc*. 2021;147:1055–1069. https://doi.org/10.1002/qj.3959

## APPENDIX A. DIFFERENCES IN THE ANALYSIS METHOD FROM T17

The following are the major similarities and differences in the datasets and analysis method of this study from those in T17.

The fundamental components of the scaling analysis method are the same in both sets of analysis, namely:

1. How a rain cluster is defined (as a set of connected rainy grids). Rain in both analyses was the total precipitation;

2. How the empirical probability distributions for A and R (as simple-scaling ansatz for a finite domain) is defined (Equation (4));

3. The hypothesis that the finite-size scaling exponent has a linear relation with the inverse of the size of the domain: $\zeta_S(z) = m_S \times \frac{1}{z} + c_S$. The thermodynamic limit of $\zeta_S(z)$, that is, $\lim_{z \to \infty} \zeta_S(z) = c_S$, is the exponent sorted.

Besides the similarity in the key assumptions of the analysis, the other details in the methods are made as similar as the differences in the two datasets would allow:

1. **Partitioning of the ocean basins (Figure 1) into domains of size Z:** We have included a figure showing the domains used in T17 (Figure A1) that can be compared with those in this work (Figure 1). Due to

the difference in the spatial resolution (36 km for WRF and 25 km for T17) of the two products, we began with $Z = 96$ for the modelled rain product. This $Z = 96$ domain covered nearly the same geographical area as the $Z = 120$ used in T17 for each ocean basin. Next we proceeded to divide this large domain into smaller and smaller “square” domains. Naturally due to the different divisibility of the largest $Z$ (96 for the current work, and 120 for T17), the possible $Z$ values than for T17 (120, 60, 40, 30, 24) are different from this work (96, 48, 32, 24). What is important to note is that the differences in the $Z$ values used between the two analyses are not really that critical; what mattered was that we had a similar range of $Z$ values in this work when compared to the T17 in estimating the exponents. Similar range of the Zs would mean that the uncertainties in the exponent estimates were about the same as those of T17.

2. **Identification of the regions which were rain predominantly from deep convection:** The idea in the original analysis was to identify regions (those square domains of the various Zs) where the rain clusters were mostly due to deep convection and not shallow convection. TRMM rainfall in T17 was the total rainfall, while we had both the grid-resolved rain from the cloud microphysics scheme and the cumulus parametrization scheme. Comparing the accepted and rejected boxes for the various Zs in Figure A1 with those in Figure 1, the boxes that were accepted for further analysis for both T17 and this work were largely similar in terms of their geographical locations. The only major difference in this work was that the method we used resulted in the rejection of all $Z = 96$ and $Z = 48$ domains, while T17 has “valid” domains for $Z = 120, 60$ and $40$ for AO. As we wish to have as close a range of $Z$ values for the linear regression for the exponents (see above), we included two additional $Z = 48$ domains (Figure 1b) which passed our test.

3. **Spanning clusters in modelled rainfall:** The discarding of the spanning clusters in the analyses were only done in the current analysis and not in T17 as there were no spanning clusters in T17. This difference was noted and discussed in the manuscript. We believed the existence of spanning clusters in model but not in observations deserved further investigation. In terms of our analysis, the removal of the spanning clusters was a necessity as we were mainly concern with the clusters that have size smaller than the domain.

4. **Empirical fitting of the probability distributions:** Both T17 and the current analysis used the maximum likelihood method to determine the exponents from the sample distributions (Equation 4).
MISR High Cloud Fraction

<figure><img src="image_18.png" /><figcaption>**FIGURE A1** Domains for T17 (their figure 2). The boxes with solid (dotted) edges are domains that passed (failed) the selection criteria imposed by T17. The contours are the satellite-derived mean annual averaged cloud fractions (reported as percentages) from 2000 to 2012 for clouds with cloud top height between 7 km and 15 km. Refer to T17 for details on the satellite-derived cloud climatology</figcaption></figure>

APPENDIX B. COMPARISON OF THE SCAL-
ING EXPONENTS FOR R AND A IN T17 WITH
THOSE IN THIS STUDY

In Figure B1 the upper panel is reproduced from Figure 4 of T17 showing the estimates for ζ<sub>A</sub> (black) and ζ<sub>R</sub> (grey)

for the four ocean basins estimated from observed rain
clusters. The lower panel is Figure 4e from this study with
the same set of estimated scaling exponents. The estimated
β values for individual oceans from this study are in Table 1
in the main text, while in T17, an estimate of β valid for all
ocean basins is 1.33 ± 0.03.

<figure><img src="image_19.png" /><figcaption>**FIGURE B1** Comparison of the scaling exponents for *R* and *A* over the four ocean basins (IO, WP, EP and AO) in T17 (upper panel) with the corresponding values in this study (lower panel)</figcaption></figure>