Estimating vertically averaged energy dissipation rate

The energy dissipation rate is an important characteristic of turbulence; however, its magnitude in observational profiles can be misidentified owing to its erratic evolution. By analysing observed data from oceanic turbulence, we show that the vertical sequences of depth-averaged energy dissipation rates have a scaling property, and propose a method to suitably estimate the vertically averaged value by utilizing that property. For scaling in the observed profiles, we found that averaging neighbouring points increases the expected value of its logarithm proportionally to the logarithm of the averaging interval. 5 Furthermore, the population mean can be estimated for the logarithm of the vertically averaged energy dissipation rate from a single observation profile, by scaling up and promoting the observed value at each depth to one that corresponds to the whole profile. The estimate allows to distinguish whether an observational profile exhibits a momentarily high value by intermittency or maintains high energy dissipation on average.


Introduction
The clear importance of determining the energy dissipation rate to study ocean general circulation has been highlighted in many studies (e.g.Munk and Wunsch, 1998;Gregg et al., 1973).However, the calculation of its average magnitude in observational profiles can be easily misled by its erratic evolution.To determine the accurate state of turbulent energy, we need to clarify the definition of energy dissipation rate.To this end, we first define statistics regarding the energy dissipation rate, which is a term conveying different meanings including the local, depth-averaged, and ensemble-averaged rate, and hence should be properly established.
By taking the expectation, the ensemble-averaged energy dissipation rate is given by where p(ǫ(x)) is the probability distribution function of ǫ(x).
If we assign a value to a vertical segment z − r 2 , z + r 2 ⊂ R in an observed profile at horizontal location x ∈ R 2 , the corresponding depth-averaged energy dissipation rate is given by ǫ r (x, ω) = 1 ǫ( x , z ′ , ω)dz ′ . (3) Energy dissipation rates obtained from oceanic observation (e.g.Itoh et al., 2009, and references therein) should be the depthaveraged ones (3) over some vertical segment with length scale r, because only the information from the segment is available, and points within the segment are indistinguishable.
Next, we discuss observational data averaging.In the study of ocean energetics, it is important to determine how turbulence energy dissipates on a large spatial scale and climate time scale, which should be represented by the 'average' energy dissipation rate (e.g., Waterhouse et al., 2014).Although previous studies have employed either arithmetic averaging, which conserves energy, or geometric averaging, which seems effective for a log-normal population (e.g., Whalen et al., 2015), no clear reasoning exists about the appropriate method of averaging observed energy dissipation rates, which would help us acquire a perspective on large-scale ocean mixing.
Geometric averaging of the vertically neighbouring values of ǫ r in Eq. ( 3) gives where Note that suffix z on the left-hand side of Eq. ( 4) is dropped because it is defined for the profile.Assuming that the data in approximately obey a Gaussian distribution for each j, ǫ r ( x j , ω ′ ) will serve as an estimate of the first parameter for the distribution of the energy dissipation rate at scale r in a profile with length ℓ = Kr.
In contrast, by arithmetically averaging the vertically neighbouring values of ǫ r in Eq. (3), we could estimate the depthaveraged energy dissipation rate for a longer vertical segment: where Observing that the ω k 's represent distinct events because we cannot simultaneously consider all K observations at z k 's, ǫ ℓ ( x j , ω ′ ) corresponds to the depth-averaged energy dissipation rate of a coarse-grained event, which contains K indistinguishable observations, instead of a single observational event.
Assuming that ln ǫ ℓ ( x ) is populated according to a Gaussian distribution at a fixed horizontal location x , we can estimate the population mean by averaging sample values ln ǫ ℓ ( x , ω m ), and the population mean is probably a key to the accurate estimation of energy balance.Hence, rather than the statistics given in Eq. ( 4) or ( 6), we should focus on the population mean for the logarithm of the vertical averages repeatedly sampled at each horizontal location.Specifically, assume that the data in approximately obey a Gaussian distribution N (µ j , σ 2 j ) for each j.The first parameter µ j for the distribution of the energy dissipation rate at scale ℓ is estimated to be which can be derived from a series of M profile observations However, we do not have access to this estimate because it requires many repeated observations at a fixed horizontal location, each one offering a vertically averaged instantaneous value ǫ ℓ ( x j , ω m ).
We propose an alternative for Eq. ( 8) to estimate the first parameter for the distribution of the vertically averaged energy dissipation rate for a profile by using a sequence of depth-averaged rates for shorter segments.To do so, we need to know the statistical relation between the depth-averaged energy dissipation rates in scales ℓ and r, for ℓ ≫ r.To bridge different scales, we unravel the scaling law that underlies the observational data.
In this regard, we use some hints from previous studies.Kolmogorov (1962) claims that the velocity correlation of homogeneous isotropic turbulence is expressed as a function of the volume-averaged energy dissipation rate.Given the intermittency of turbulence, the logarithm of the depth-averaged energy dissipation rate obeys a Gaussian distribution, whose expected value and variance depend linearly on the logarithm of the averaging scale.This claim suggests that the profiles of depth-averaged energy dissipation rates taken from any turbulence can have scale dependencies in general.
In this paper, we show that the observed data of depth-averaged energy dissipation rate have a scaling property and formalize a method to estimate the first parameter for the distribution of the vertical average by using this property.The remaining of the paper is organized as follows.In the next section, we analyse the scaling property and propose a calculation method for the average.Then, we present calculation results by using real observational data.In the discussion, we describe the interpretation and importance of these results, and close the paper with conclusions on the importance of the proposed method and perspective for applications.

Proposed method
In this section, we describe the turbulence observation data employed in this study, formalize the scaling dependence encoded in the data, and detail a method for calculating the vertical average accounting for the scaling relation.

Observational data
The turbulence observational data were retrieved from the Pacific ocean (Fig. 1) (Goto et al., 2018) and comprise approximately 400 profiles, each of which typically extends over 2000 to 6000 m-depth below the sea surface, with observational points every 5 to 10 m.
These positive-valued data have the following characteristics: 1.For sequential data in ordered set each profile exhibits a very erratic evolution that impedes the recognition of a continuous curve along the depth direction (Fig. 2(a)).
2. After taking the logarithm of the values, the sequences seem more continuous (Fig. 2(b)).If we treat individually all the available data in nonordered set the distribution of the logarithmic values approximately resembles a Gaussian distribution (Fig. 4), as pointed out by Gregg et al. (1973).We use this log-normality only in the approximate sense.
3. For the logarithmic sequence in ordered set after averaging each vertically neighbouring pair of values ǫ r in a profile (Fig. 2(c)), the sequence shows larger values of ln ǫ r ′ than the original sequence on average: where E z represents averaging in any vertical interval larger than r ′ (Fig. 2(d)).
4. If we repeat procedure 3, we obtain the sequences that have even larger values of ln ǫ r ′ on average (Fig. 3).
The above shift is a direct consequence of the inequality of arithmetic and geometric means, i.e.
for the upper, (U ), and lower, (L), parts of ǫ 2r .However, the shift amount can be a specific property of the analysed data.
Furthermore, although the shifts ln ǫ 2r ′ − ln ǫ r ′ are fluctuating from point to point, it is suggested that the vertical average of the shifts for a profile does not differ much between all the observed profiles.

Scaling analysis
The statistical phenomenon observed in Fig. 3 can be attributed to the scaling properties encoded in the data sequence from the profiles.Here, we introduce a random variable α, which serves as the factor of proportionality, and postulate that averaging neighbouring points increases its logarithmic value proportionally to the logarithm of the averaging interval: where the right-hand side is αn ln 2 for r ′ = 2 n r.By taking the expectation, we obtain The target population here is all possible states of profiles at any horizontal points of observation.Thus, if we interpret E as averaging E x among all the available pairs of segments with length r and r ′ taken from observed profiles at various horizontal locations, then Eq. ( 10) serves as an estimator for the factor: Note that the estimated value α is accompanied by an estimation error, which will be discussed later in section 3.1.2.
If this assumption holds and factor α is calculated from observational data, then we obtain an estimator for the vertically averaged energy dissipation rate in a larger segment: which corresponds to several ln ǫ r ( x , z)'s centred at z satisfying data satisfy where the geometric summation is taken for every pair of segments with lengths r and r ′ taken from profiles at various horizontal locations.
If the observational points in a profile are arranged with uneven spacing and their coordinates are approximated by a power function z k = ak c + b, k = 0, 1, 2, • • • , the left-hand side of Eq. ( 9) should be corrected as c −1 ln (ǫ r ′ /ǫ r ) for the right-hand side to share the common form αn ln 2. In detail, the scaling analysis based on Eq. ( 9) is performed as follows: 1. Suppose we have an observational sequence in ordered set for each j.
2. By averaging vertically neighbouring 2 n points, we obtain ordered set where r ′ = 2 n r and As a typical observational segment is r ≃ 5 m, with profile length 5000 m, n is not larger than 10 ln (5000/5)/ ln 2. 3.With respect to all possible triplets (j, k, k ′ ) for fixed n, we calculate the average of where c j is the correction for the uneven observational spacing, and denote the average as Y n .

Rescaled average
Suppose S( x ) represents the observation of a profile at coordinate x , i.e. the set of all pairs of segments with length ℓ and r in a profile of length ℓ = Kr.Assuming statistical homogeneity of α in Eq. ( 9) and applying it to the events in S( x ), the conditional expectation is This implies an estimator for the first parameter for the distribution of vertical average ǫ ℓ ( x ): where we used ln ǫ ℓ ( x , z ′ ) = ln ǫ r ( x , z) + α ln K for the third equality.Here, 'the first parameter' is for the log-normal distribution of the states that can occur during the profile observation.Note that in general, estimate ǫ ℓ ( x ) is neither equal to geometric mean ǫ r ( x ) in Eq. ( 4) nor to arithmetic mean ǫ ℓ ( x ) in Eq. ( 6).Instead, it is similar to ǫ ℓ ( x ) in Eq. ( 8) except that every ln ǫ ℓ ( x , ω m ) is substituted by the estimated ln ǫ ℓ ( x , ω k ), which resemble K repeated observations (see Fig. 5).Hereafter, the estimate given by Eq. ( 15) is called the rescaled average.
We cannot necessarily assume statistical homogeneity of α along the horizontal or vertical direction.With this in mind, we also consider the conditional expectations using three different α values: where i denotes either three horizontal regions, as those divided by dashed lines in Fig. 1, or two vertical regions above and below 200 m depth.Consequently, α i is derived from the scaling analysis using the data from region i, and thus we have other approximation formulae for ǫ ℓ ( x ): according to the region that x belongs to. 3 Results

Estimated value
By using the observational data, we conducted the scaling analysis according to the procedure described at the end of section 2.2.The values of E x [ln (ǫ r ′ /ǫ r )] are plotted against ln (r ′ /r) as black dots in Fig. 6, where a scaling relation with expo-5 nent α = 0.467 valid from r = 5 m to r = 5000 m can be observed.
The comparison among data above and below 200 m shows that the scaling appears to be statistically homogeneous along the depth direction (top graph of Fig. 6).On the other hand, the comparison among data from the three horizontal areas does These results confirm the claim that averaging neighbouring points increases the expected value of its logarithm by a rate, α > 0 or α i > 0, proportional to the logarithm of the averaging interval.Moreover, we determined the values of the scaling exponent through the statistical analysis of the observational data.

Estimation error
For simplicity, here we only consider the case with a single α.The scaling relation ( 11) provides a mean estimate for α from the data as We can also calculate the variance of the estimated α from the data as where σ 2 intra is the average of variances within a profile, and σ 2 inter is the variance of profile-averages.Figure 7 shows that σ intra and σ inter are approximately 0.44 and 0.20, respectively, when r ′ ≫ r.

Estimated value
The vertical average of the energy dissipation rate for the upper ℓ = 2000 m at each horizontal location was calculated according to Eq. ( 15) or ( 17).The value of ℓ was selected such that the length of the calculated segments agrees with that of a typical observed profile.We excluded the profiles with less than 50 data points in the upper 2000 m.
In Fig. 8, the simple arithmetic average from Eq. ( 6) is compared with the rescaled estimates using a single α value (top graph) and three α values (bottom graph).In most points, the rescaled estimate shows larger values than the arithmetic average, which still exhibits very large values intermittently, unlike the rescaled estimate.There seems to be no significant difference between the distributions in the top and bottom graphs of Fig. 8.These results suggest that we can suitably estimate the order of the values even using a single α value.inter is the variance of profile-averages.

Estimation error
Here, we only consider the case with a single α for simplicity, and fix a horizontal location x .We assume that ln ǫ ℓ ( x ) is populated according to an almost Gaussian distribution with unknown mean µ and variance σ 2 , which are to be estimated.
Using the arithmetic averaging strategy, we obtain a single sample ln ǫ ℓ = ln K −1 K k=1 ǫ r ( x , z k ) with variance σ 2 ln ǫ ℓ = σ 2 , and the likelihood of population mean µ is proportional to exp (µ − ln ǫ ℓ ) 2 /(2σ 2 ln ǫ ℓ ) .Note that, basically, we cannot obtain any information about the estimation error.
By contrast, using the rescaled averaging strategy, we retrieve K samples Inserting typical values in the observed data for K ≫ 1, we obtain where σ 2 ln ǫr = 2.6 2 is the average of intra-profile variance of ln ǫ r calculated from observational data, and σ intra and σ inter are from sect 3.1.2.In particular, for K = 2000 m/5 m = 400, we have σ ln ǫ ℓ = 1.2, which is 0.52 in common logarithm.
Comparing the estimation errors σ ln ǫ ℓ and σ ln ǫ ℓ , although the rescaled average is not necessarily more accurate than the arithmetic average, the former has explicitly defined uncertainty that remains relatively unchanged from profile to profile, which is an advantage as an estimate for the first parameter.

Horizontal distribution
Figures 9 and 11 show the horizontal distribution of the vertically averaged energy dissipation rate by using a single and three α values, respectively.A cross section along 47 • N is shown in Fig. 10.The distributions of the areas with large ǫ ℓ value are different in many regions between the rescaled and arithmetic averages, despite some similarities.High peaks are sometimes isolated for the arithmetic average, whereas they tend to form clusters for the rescaled average (see Fig. 10), which has at most the order of 3 × 10 −8 to 1 × 10 −7 W • kg −1 (see Figs. 9 and 11), suggesting that the rescaled average captures more robust features from the distribution of energy dissipation.
These results indicate that the calculated values from rescaled averaging can be used to correctly interpret the horizontal distribution of energy dissipation.

Discussion
1. Previous studies have not clarified the relation between small-scale observation of turbulence in several meters with data substantially intermittent and the large-scale features of energy dissipation, which should be more persistent.In this paper, we propose a method to relate them quantitatively by elucidating the scaling relations underlying the observational data and thereby estimating the vertical average through rescaled averaging.
2. The volume-averaged energy dissipation rate can be considered as the value when the local energy dissipation rate, a Schwartz distribution, is evaluated against a test function φ: where the integration is a formal notation.For example, expected value E ln ǫ 2 n r (z0) can be regarded as the comparison between the estimation of local energy dissipation rate ǫ(z 0 + z) and its rescaled version where we use indicator function φ(z) = r −1 ½ [−r/2,r/2] (z) as test function.In terms of the local energy dissipation rate, the scaling law in Eq. ( 9) means because 3. In principle, the scaling property represents the self-similarity underlying the local energy dissipation rate defined along the depth axis.
From Eq. ( 10), E ln 2 −αn ǫ 2 n r (z k ) ǫr(z k ) = 0, which implies that the sequence of 2 −αn ǫ 2 n r (z k ) is similar to that of ǫ r (z k ) in the logarithmic scale.We can illustrate the self-similarity inherent to the profile data by scaling them appropriately to obtain similar sequences ζ 2 n r in Fig. 12. Specifically, if we use a test function defined with the aid of indicator function and rescale it to the evaluation of ǫ against this test function is given by which has scale 2 αn .We can thus define statistically invariant sequence around z 0 for each n = 0, 1, • • • .

Conclusions
1.In this paper, we show that the logarithmic sequences of observed depth-averaged energy dissipation rate exhibit a scaling dependence, and formalize the estimation of the population mean for the logarithm of the vertically averaged rate by using this dependence.
2. From the scaling property in the observed profiles, we found that averaging neighbouring points increases the expected value of its logarithm by a rate proportional to the logarithm of the averaging interval.Furthermore, this property allows us to estimate the population mean for the logarithm of the vertically averaged energy dissipation rate from ocean turbulence profile data.
3. By scaling up and promoting the observation value at each depth to one corresponding to the whole profile, information from a virtual series of profile observations can be extrapolated from a single profile observation.
4. The estimate of the vertically averaged energy dissipation rate is accompanied by an explicitly defined estimation error for the population mean, and serves as a measure for distinguishing momentarily large but intermittent values from persistently large energy dissipation on average.In addition, the estimated vertical average should correspond to the energy dissipation rate used in ocean general circulation models, in which case the water column depth can be considered for ℓ in Eq. ( 16).
5. The scaling exponent depends on the horizontal location, as the turbulence differs among places.Therefore, more observations are needed to establish a more accurate scaling relation.To estimate the spatial distribution of energy dissipation, data assimilation into an ocean general circulation model can be performed using the procedure proposed in this study.

Appendix A: Measurement error in scaling analysis
If observational data ǫ obs,r contain measurement error η r obeying N (0, Σ 2 ), from logarithmic manipulation we obtain This has the following effect on scaling: where D r ≡ 1 2 Σ 2 E (ǫ r ) −2 and h(−2) is the scaling exponent for ǫ −2 r .The data may contain errors at most in the order of D r > 0. Hence, we use a curve of the form y = ax + b(1 − e −cx ) for fitting, and regard c as the slope for ln (ǫ r ′ /ǫ r ).

Nonlin.Figure 1 .
Figure 1.Horizontal locations of the observed profiles (red), and land-sea boundaries (green).The units of longitude and latitude are • E and • N, respectively.

Figure 3 .
Figure 3. Evolution of averages Ez [ln ǫ r ′ − ln ǫr] as r ′ increases.For simplicity of presentation, Bèzier curves are used for average Ez.

Figure 4 .
Figure 4. Distribution of logarithm of observational data (bar chart), and the best-fitting Gaussian distribution (red).

Figure 5 .
Figure 5. Concept of rescaled average.A single observation (green square) is rescaled to the profile length (purple rectangle) and treated as an ensemble member.

Figure 6 .
Figure 6.Scaling property of observational data.Average of ln (ǫ r ′ /ǫr) (dots) and best-fitting curve of the form y = ax + b(1 − e −cx ) using all the data points (solid line).The results were obtained for different depth sections (top graph), and horizontal sections (bottom graph).

Figure 7 .
Figure 7. Estimation errors of α according to scales ln (r ′ /r).Intra-profile σ 2 intra is the average of variances within a profile, and inter-profile σ 2inter is the variance of profile-averages.

Figure 8 .Figure 9 .Figure 10 .Figure 11 .Figure 12 .
Figure 8. Scatter plot of depth-averaged energy dissipation rates for upper 2000 m comparing the expected value obtained from rescaled and arithmetic averaging.The top and bottom graphs show the results using a single and three α values, respectively.19