Reporting delays: A widely neglected impact factor in COVID-19 forecasts

Abstract Epidemic forecasts are only as good as the accuracy of epidemic measurements. Is epidemic data, particularly COVID-19 epidemic data, clean, and devoid of noise? The complexity and variability inherent in data collection and reporting suggest otherwise. While we cannot evaluate the integrity of the COVID-19 epidemic data in a holistic fashion, we can assess the data for the presence of reporting delays. In our work, through the analysis of the first COVID-19 wave, we find substantial reporting delays in the published epidemic data. Motivated by the desire to enhance epidemic forecasts, we develop a statistical framework to detect, uncover, and remove reporting delays in the infectious, recovered, and deceased epidemic time series. Using our framework, we expose and analyze reporting delays in eight regions significantly affected by the first COVID-19 wave. Further, we demonstrate that removing reporting delays from epidemic data by using our statistical framework may decrease the error in epidemic forecasts. While our statistical framework can be used in combination with any epidemic forecast method that intakes infectious, recovered, and deceased data, to make a basic assessment, we employed the classical SIRD epidemic model. Our results indicate that the removal of reporting delays from the epidemic data may decrease the forecast error by up to 50%. We anticipate that our framework will be indispensable in the analysis of novel COVID-19 strains and other existing or novel infectious diseases.

The COVID-19 pandemic has been crippling the world's health, economies, and quality of life for over three years.While we are likely beyond the peak, we are gradually understanding that COVID-19 is here to stay.Fast mutation rates of the virus and its overwhelming spreading capability make it extremely hard if not impossible to eradicate [1].COVID-19 is not the first and almost surely not the last pandemic to hit humanity.Therefore, in order to better prepare and withstand other contagious diseases in the future, we need to extract as many lessons as possible from the COVID-19 pandemic.
Public awareness is, arguably, the first line of defense against any infectious diseases.Efficient collection of epidemic data and accurate epidemic forecasts allow for timely containment of the spread or flattening of the curve to win time for the development of pharmaceutical treatment methods.Due to the success of network epidemiology and the broad availability of data, we have made significant advances in epidemic modeling and forecast methods [2][3][4][5][6].However, the accuracy of epidemic forecasts strongly depends on the accuracy and timeliness of the input data.
Are epidemic data -especially in the short-time period after the onset of the epidemic -devoid of noise?Multiple sources suggest the negative answer.One problem is delays in data reporting [7,8], defined as the time difference between the event when the person was affected by the virus and the time this event is accounted for.The reporting delays may occur due to many reasons, including patient hesitancy, medical testing delays, and reporting delays by public health authorities [9,10].As a result, reporting delays may vary not only across different diseases but also across different regions of interest.The median diagnosis delay for Malaria is, for instance, approximately four days [11].Similarly, the research on the Middle East Respiratory Syndrome CoronaVirus (MERS-CoV) found a time difference between the symptom onset and confirmation of approximately four days [12].Hepatitis A, Measles, and Mumps data are usually reported after eight days [9].Reporting delays for some diseases are significantly longer: Hepatitis B, Shigellosis, and Salmonella can take on average 2-3 weeks [9].Recent COVID-19 studies reveal significant reporting delays of infections in China [13][14][15][16][17][18], Italy [19], Germany [20,21], Singapore [22], the USA [23], and the UK [7].
In our work, we develop a statistical framework to de-noise epidemic data by removing reporting delays.We demonstrate that the removal of reporting delays may significantly improve the accuracy of epidemic forecasts.Our statistical framework can be used in combination with any epidemic forecast method that takes infectious, recovered, and deceased epidemic data as input.Our work is organized as follows.We first present the evidence for the presence of reporting delays in the epidemic data extracted from the first COVID-19 wave.We then proceed to develop a statistical framework to remove reporting delays.After validating our framework on synthetic data, we move on to remove and analyze reporting delays in 8 hotspots of the COVID-19 pandemic.We conclude our work with a discussion of the impact of delay removal on the accuracy of epidemic forecasts.

A. Notation
Before presenting our findings, we introduce our notation for the epidemic data.Throughout the text, we operate with the infectious I, recovered R, and deceased D data.Each dataset is a time series of values, each corresponding to a specific observation date.For brevity, we refer to the triplet of infectious, recovered, and deceased data as Y = {I, R, D}.All values contained in the Y time series are fractions of individuals found in the corresponding state on a specific day.For instance, I[k] corresponds to the ratio of individuals who are infectious on day k and the population size.Since COVID-19 reported data often comes in the form of changes in the number of epidemic cases, we find it convenient to introduce the daily changes in epidemic data as Further, in this work, we operate with reported epidemic data Ỹ and inferred data Ŷ .Since reported and inferred data are expected to differ from the true data, we need to distinguish the three.We refer to the true data as Y .We summarize our notation in Table I.

B. Evidence for reporting delays in epidemiological data
We begin the exposition by considering daily reports on the number of infected Ĩ, recovered R, and deceased D individuals in Spain.While both the infected and the deceased data indicate that the first COVID-19 wave in Spain peaked in April 2020, the exact timings of the two peaks are, nevertheless, different.As observed in Fig. 1a, reported new deceased cases ∆ D reached their peak on April 1st, 2020, while the highest fraction of infectious individuals Ĩ was observed 22 days later on April 23, 2020.This observation is not specific to Spain: the peaks in the number of reported new deceased cases precede those of infectious cases by more than one week in most regions, Fig. S1.We make similar observations for COVID-19 daily recovery reports ∆ R, which exhibit their peaks after those of the deceased cases ∆ D, Fig. 1a, and Fig. S1.Furthermore, when plotted as a function of daily deceased cases ∆ D, infectious cases Ĩ and daily recovered cases ∆ R form loop patterns, see Fig. 1b,c and Fig. S2.The patterns observed in Fig. 1a,b,c, Fig. S1, and Fig. S2 may indicate the presence of reporting delays.Indeed, the SARS-CoV-2 virus has hardly changed during the first wave of the pandemic.Therefore, recovery γ r and death γ d rates are estimated to be approximately constant during this period [24][25][26], implying that changes in the fractions of recovered ∆ R, and deceased ∆ D data are proportional to the fraction of infectious individuals I, contradicting Fig. 1a,b,c.We hypothesize that the disagreement between ∆ R, ∆ D, and I is due to reporting delays.If each ∆ R, ∆ D, and I time series are reported with different delays, their peaks are expected to differ.Likewise, the loop patterns of Fig. 1b,c may also be the result of an effective time shift between two non-monotonous time series.The upper part of the loop in Fig. 1b is due to the fact that initially ∆D[k] and I[k] both increase as a function of time.The lower part of the loop is observed when ∆D[k] and I[k] are both decreasing as a function of time.The middle section of the loop pattern corresponds to the time window when ∆D[k] is decreasing while I[k] is still increasing as a function of time.We verified our hypothesis on synthetic epidemic data with reporting delays, observing similar loop patterns, see Supplementary Information B and Fig. S3.
C. A null model for reporting delays reporting delay reported time FIG.2: A schematic representation of epidemic events and corresponding delay times.
To uncover reporting delays, we employ the following null model.Within the model, each individual i is endowed with random event time τ Yi of either getting infected, τ Ii , recovered, τ Ri , or deceased, τ Di .We denote the random delay time for event and the corresponding random reported time as Assuming that the reporting delays T Yi are independent of events Y i , we obtain for the reported time By taking the arithmetic mean of Eq. [1], we find that the expected reported fraction of individuals in state where Y = {I, R, D}.Equation [2] serves as the foundation for uncovering reporting delays and improving epidemic forecasts.

D. Statistical framework to uncover reporting delays
While the delay inference problem can be formulated and solved in a non-parametric way, to simplify the exposition, we solve the problem parametrically.We assume that the shapes of delay distributions Pr[T Y = m] are known, but their parameters Y are not.We then search for κ Y maximizing the correlations among the epidemic data: ∆Y = {I, ∆R, ∆D}.
We consider three families of two-parameter discrete-time distributions: the negative binomial distribution, the Neyman type A distribution, and the Pólya-Aeppli distribution, see Materials and Methods A. These distributions are quite general and contain some well-known 1-parameter distributions as special cases.The logarithmic, geometric, and Poisson distributions are, for instance, all special cases of the negative binomial distribution [27].The three distributions are sufficiently general to generate data with variable mean and variance and of varying skewness [27,28] -properties expected of the epidemic delay data, see Materials and Methods A. For brevity, we focus on the Pólya-Aeppli distribution in the main text and report the results for the other two distributions in Fig. S4.
In more precise terms, we assume that the reporting delays correspond to three datasets, Y = {I, R, D}, which are all characterized by the Pólya-Aeppli distribution, albeit with different parameters κ Y ≡ {λ Y , θ Y }.Hence, there are in total 6 parameters for three Pólya-Aeppli distributions, κ = (λ For prescribed values κ, we can remove the effects of reporting delays by solving Eq. [2] as a system of linear equations for original epidemic data ∆Y .Since the resulting epidemic data is expected to differ from the true data, we refer to the reconstructed data as ∆ Ȳκ .Given the incremental time series ∆ Ȳκ , the cumulative time series Ȳκ follows as In our framework, we aim to determine parameters κ that maximize the product of pairwise correlations among the three epidemic time series in Y : where ρ(X, Y ) is the Pearson correlation coefficient between time series X and Y .The objective function O b (Y ) in Eq. [3] reaches its maximum value of 1 when all three pairwise correlations among the I, ∆R, and ∆D time series are 1, which we expect when recovery γ r and deceased γ d epidemic rates are constant.Due to the nature of the Pearson correlation coefficient, As a result, we can only infer the reporting delays up to a constant time shift T .
To maximize O b (Y ), we use the random search [29]: we conduct a large set of independent random iterations, = 1, . . ., L. At each iteration , we select the elements of the parameter vector κ uniformly at random from the prescribed domain of values.Further, at each iteration , we use the selected parameter set κ to reconstruct the original epidemic data Ȳκ by solving Eq. [2].We then compute the pairwise correlations in Ȳκ to obtain the objective function O b Ȳκ .After completing all random search iterations, the resulting delay parameter κ is the one maximizing the objective function, κ = argmax κ O b Ȳκ , see Materials and Methods B.
where β, γ r , and γ d are model parameters quantifying, respectively, the infection rate, the recovery rate, and the deceased rate.
To test our framework, we generated 50 SIRD epidemic model datasets with different epidemic parameters and added synthetic reported delays to the obtained times series, as prescribed by Eq. [2].After generating the synthetic data we forgot the parameters of reporting delays and used our statistical framework to infer them.Figures 3a,b indicate that the inferred delay parameters are in good agreement with the true parameters that generated the datasets.We observe that the inference errors decrease fast as a function of the number of iteration steps saturating at the margin of 2 days after 10 7 iterations, Fig. 3a.To assess the robustness of the inference procedure, we have conducted cross-inference experiments by generating synthetic data with one delay distribution and inferring delays using another distribution, arriving at similar results, Fig. S4.

F. Uncovering reporting delays in real data
After testing our inference framework on synthetic data, we moved on to uncover reporting delays in real epidemic data that we have collected from 8 regions worldwide, Supplementary Information A. Figs.1d-f display the reconstructed epidemic data for Spain after the identification and removal of reporting delays, see Table S2 for the inferred parameters.
Consistent with our expectations, the reconstructed infected, recovered, and deceased time series are strongly correlated, Figs.1e,f, and their peaks co-occur at day 40, Fig. 1d.By comparing the original correlation coefficients O b ( Y ) with their optimized counterparts O b ( Ȳ ) we discover that the latter are significantly increased with a relative improvement ranging from 0.94 to 9.5 after the reconstruction, Table S2.
We summarize the properties of uncovered reporting delays in the 8 regions in Table 2.We find that the infectious and recovered data delays are longer than those for deceased cases and vary from several days to several weeks.This observation is hardly surprising.Indeed, there are several factors contributing to the delays of infectious cases.One factor is the delay between an individual becoming infectious, and the symptom onset [30].Another factor, particularly significant in the first COVID-19 wave, is the delay between the symptom onset and the test result [31].In their turn, delays in the recovered events are likely caused by hospital discharge policies.This is the case since the recovered data is usually derived from hospital discharge events, which occur after the patient's recovery.
Based on the expected delay values, one can naturally split the 8 ROIs into three categories: (i) large infectious and small recovered delays: Romania, Germany, and Denmark, (ii) small infectious and large recovered delays: Italy and Spain, and (iii) large infectious and recovered delays: Wuhan, Hubei, and Turkey.
Small infectious delays imply that COVID-19 tests are timely and accurate.This seems to be the case for Italy, which executed more tests per capita in April 2020 than other countries [32].On the contrary, the testing ability of the Hubei province was significantly insufficient during the first wave of the COVID-19 outbreak.
Large delays in the recovered data may be attributed to strict hospital discharge policies.Discharge policy in Italy, for instance, was based on the negative test [33], and it has been shown that COVID-19 tests may stay positive for an extended time after COVID-19 symptoms disappear.In contrast, discharge policies in Denmark were based not on the negative test but on patient symptoms [33], likely leading to shorter reporting delays in the recovered data.
Large standard deviations in the reporting delays may indicate irregularities in reporting mechanisms.As an example, Hubei province expanded its daily testing capacity from 200 to 2,000 individuals from the beginning of the pandemic until January 27th [34].As a result, fewer individuals were tested late at the end of the first COVID-19 wave compared to its beginning, likely resulting in the large standard deviation of infectious data delays.
The small standard deviation in the delays of recovered data observed in Hubei and its capital Wuhan is likely the consequence of the strict discharge criteria [35].Although strict discharge policies cause significant delays, these delays are similar, resulting in relatively smaller σ R values.In contrast to Chinese regions, the recovery data for Germany are not reported directly but instead are estimated by a not explained algorithm [36], resulting in large errors and, consequently, larger σ R values.Based on the optimized value O b ( Ȳ ) as shown in Table 3, the optimization performances of our algorithm for Italy, Spain, and Turkey are better than the other countries.Fraction of recovered data Uniformly at random select delay parameters.

G. Improving epidemic forecasts
Is accounting for reporting delays likely to improve epidemic forecasts?To answer this question, we designed two experiments.Experiment 1 aims to forecast the epidemics ignoring reporting delays and serves as a baseline for Experiment 2, which uncovers reporting delays prior to forecasting epidemic data.
In both experiments, we split the reported epidemic data into two parts, which we call the training and the testing sets, respectively, see Fig. 4a.In Experiment 1, we fit the training set with the SIRD epidemic model, obtaining model parameters β, γ r + γ d , and the fraction of initial infected cases I[0].We then use the SIRD model with the obtained parameters to forecast the epidemic data, Supplementary Information C.
Regions  S1 and S2 for the corresponding parameters of the Pólya-Aeppli distributions.
In Experiment 2, we first use the reported data Y in the training set to infer the parameters of the delay distributions κ.We rely on these parameters to remove reporting delays from the testing set and obtain reconstructed data ∆ Ī.In the next step, we fit the reconstructed data ∆ Ī with the SIRD model, obtaining β, γ r + γ d and I[0] spreading parameters.We use these spreading parameters to forecast the epidemic data ∆ Î, which we compare to the reported data ∆ I in the testing set.Since ∆ I in the testing set contains the reporting delays, while the forecast ∆ Î does not, we added the inferred reporting delays back to the ∆ Î using Eq.[2], obtaining ∆ Î * , See Fig. 4a and Supplementary Information C. Figures 4b,c present the results of the forecast experiments for Spain, indicating that correcting for reporting delays does improve the forecast accuracy.To evaluate the forecast errors in a systematic way, we evaluated the root mean square errors between the forecasts I F and the testing I T sets: where n is the size of the testing set.Further, to quantify the benefits of accounting for reporting delays, we used the ratio of the root mean square errors measured in forecasts with and without reporting delays, ≡ , where Î1 and Î * 2 are the epidemic forecasts obtained in experiments 1 and 2, respectively.The smaller the ratio, the smaller the relative forecast error.
We measured the accuracy of epidemic forecasts with different start dates for Spain in Germany.While we did observe a substantial improvement in forecast accuracy for Spain, Fig. 4d, this was not the case for Germany, where the removal of reporting delays only improved the forecasts by a small margin, Fig. 4e.On a broader scale, we observed that the benefits of correction for reporting delays vary across all regions of interest, Fig. 4f, Fig. S5, and Fig. S6.While there is nearly a two-fold improvement in the forecast accuracy for Denmark, there is little improvement for Wuhan and Hubei.
There are two factors that may hinder epidemic forecast accuracy.The first factor is the insufficient accuracy of reporting delay removal.The second one is the inability of the SIRD model to reproduce the COVID-19 dynamics accurately.We can quantify the former by the maximum attained value of the objective function O b ( Ȳ ) max that we aim to maximize when removing reporting delays.While there is no direct way to quantify the goodness of the SIRD model in epidemic forecasts, as an indirect measure we use the SIRD fitting error obtained by fitting the training set with the SIRD model.As seen in Fig. 4g, the forecast error is strongly correlated with the SIRD fitting error, Pearson r = 0.89, p < 10 −27 , indicating that the highest epidemic forecast accuracy is attained when the model fit is accurate.When it comes to the benefits of removing reporting delays, we observe that the forecast error ratio depends both on the O b ( Ȳ ) max and SIRD model fitting, Fig. 4h.We observe the largest error ratio in Wuhan, Hubei, and Romania, Fig. 4h.Wuhan and Hubei correspond to the largest SIRD fitting errors, while Romania corresponds to the lowest SIRD fitting errors.At the same time, all three regions are quantified by the lowest O b ( Ȳ ) max values.The other 5 regions, corresponding to smaller error ratios, are characterized by significantly larger O b ( Ȳ ) max values, Fig. 4g and Table 2.
These observations indicate that O b ( Ȳ ) max may be used as an early indicator of the benefit of delay removals.Indeed, lower O b ( Ȳ ) max values may indicate that reporting delays were not removed successfully or the initial assumption of the proportionality among I, ∆R, and ∆D times series does not hold.We note that O b ( Ȳ ) max is a good indicator regardless of the SIRD fitting error.In the case of Romania, SIRD fitting errors are among the smallest resulting in low forecast errors.Yet, the removal of reporting delays does not result in even lower forecast errors.
There might be multiple reasons for suboptimal reporting delay removals in the case of Hubei, Wuhan, and Romania.One possibility is the non-stationary nature of delays.The main assumption of our framework is that delay times are independent and drawn from the same distribution.While COVID-19 did not significantly mutate over the short time span of the first wave, our ability to handle the infections improved significantly.The PCR test capacity in China has grown remarkably during the first COVID wave, possibly explaining the limited effect of reporting delays on improving epidemic forecasts in Hubei and Wuhan.

II. DISCUSSION
Data delays are ubiquitous in data sciences and adversely affect data analysis [37].The unique property of the epidemic data, enabling us to identify and filter out reporting delays, is the proportionality between the number of infectious individuals I and the rates of change in the deceased ∆D and the recovered ∆R individuals.
We relied on this proportionality property to develop a parametric statistical framework to uncover reporting delays in the first COVID-19 wave and applied it to 8 regions significantly affected during the first COVID-19 wave.The character of uncovered delays varies across studied regions and can be explained by region-specific medical capacities and epidemic policies.
Concerning the epidemic forecasts, we found that the benefits of curated epidemic data, as opposed to the raw data, are maximized in situations when reporting delays are removed efficiently.One of the main factors hindering the efficiency of delay removal -the non-stationarity of reporting delays -can be caused by either varying spreading properties of the virus or by rapid changes in medical capacities or epidemic restrictions across the regions.While the spreading properties of COVID-19 did not change significantly during the first wave of the pandemic, both medical capacities and epidemic policies were the subjects of constant updates as the society learned how to handle the pandemic.
In conclusion, our statistical framework does not make strong assumptions about the spreading mechanisms of a pathogen.Be it a well-studied viral infection or a novel virus, our framework should be instrumental in uncovering reporting delays, as long as the spreading properties of the virus remain unchanged during the observation period.Our framework can be used as a preliminary filter for any epidemic forecast tool that takes infected, recovered, and deceased data as input.Accurate and epidemic forecasts are of immense value for society and policymakers to minimize the adverse effects of the virus.

A. Types of distributions
To determine the family of reporting delay distributions that best suit our data, we consider three different two-parameter discrete distributions [38] The Neyman type A distribution with parameters ξ > 0 and µ > 0 has mean value E[T ] = ξµ and variance V ar[T ] = ξµ(1+µ).

B. Inferring reporting delays
To uncover reporting delays in epidemic data, we determine parameters of delays distributions κ Y maximizing the pairwise correlations between the epidemic time series, Eq. [3].This optimization problem can be compactly written as: Here T is the size of epidemic time series, Y ≡ { I, R, D} and Y ≡ {I, R, D} are reported and reconstructed epidemic data, respectively, while ∆Y = Ψ −1 Y ∆ Y are the matrix solutions of Eqs.[2].Indeed, Eqs. [2] can be written in the matrix form as ∆ Y = Ψ Y ∆Y for Y = {I, R, D}, where Pr[T Y = i − j] if i j. 0 otherwise., Then Ψ −1 is the inverse of Ψ Y and ∆Y = Ψ −1 Y ∆ Y .In the main text, we assume that reporting delays are iid random variables drawn from three Pólya-Aeppli distributions, Eq. [7], with distinct parameters {λ I , θ I }, {λ R , θ R }, and {λ D , θ D }.In the case of Pólya-Aeppli distributions, the parameter vector takes the form of κ ≡ {λ I , θ I , λ R , θ R , λ D , θ D }.We solve the optimization problem given by Eq. [9] and the random search [29].For each κ we use Eq.[2] to reconstruct original data Y , which we then use to compute the objective function O b (Y ).After completing all iterations, the thought parameter vector κ describing reporting delays is the one corresponding to the maximum O b (Y ) value, κ = arg max κ O b (Y ).

FIG. 1 :
FIG. 1: Evidence for reporting delays in COVID-19 epidemic data.a, Reported infectious Ĩ[k], recovered ∆ R[k + 1] and deceased ∆ D[k + 1] cases for the first COVID-19 wave in Spain.The k = 0 day corresponds to February 25, 2020.Note the difference in peaks of the reported data.Panels b and c display pairwise color-coded scatter plots of ∆ D[k + 1] vs Ĩ[k], and ∆ R[k] vs ∆ D[k] for Spain.Colors, from light to dark green, reflect different days in the data ranging, respectively, from k = 0 to k = 80.Note that the scatter plots in panels b and c form, respectively, clockwise and counterclockwise loop patterns.Panels d-f display the epidemic data in Spain after the removal of reporting delays.
E. Tests on synthetic dataBefore uncovering reporting delays in real epidemic data, we test our framework on synthetic datasets, which we generate with the SIRD compartmental model.Within the SIRD model, the population is split into four compartments: susceptible S, infectious I, recovered R, and deceased D. Compartment S denotes the fraction of susceptible individuals who can be infected by infectious individuals.Compartment I denotes the fraction of individuals who have been infected but have not recovered or are deceased.Compartments R and D are respectively the fractions of individuals who have recovered or are deceased.The SIRD model assumes that recovered individuals become immune and cannot be infected by the virus in the future.Discrete-time transitions between the compartments are governed by first-order difference time equations data Y = {I, R, D}.

FIG. 3 :
FIG. 3:Uncovering reporting delays in synthetic data.We generate epidemic data using the SIRD model, Eq.[4] with parameters β = 0.23, γr = 0.079, and γ d = 0.11 and variable delay parameters.The schematics on top illustrates the process of generation and inference on synthetic data.a, The mean absolute inference errors as functions of the number of iteration steps in the random search.Due to the nature of the Pearson correlation coefficient, it is only possible to infer relative delays, E [TR] − E [TD], E [TI ] − E [TD], and E [TR] − E [TI ]. b, The relationship between the true and the inferred values of delay parameters.(c,d) An example of an SIRD synthetic data (c) before and (d) after the removal of synthetic reporting delays.

FIG. 4 :
FIG. 4: Accounting for reporting delays improves epidemic forecasts a, The schematic diagram for the two forecast experiments.b, c display the results of epidemic forecasts in Spain with different forecast start dates.The forecasts of experiment 2 (yellow) are closer to the reported data (red) than the forecasts of experiment 1(blue).The forecast results for all 8 countries or regions are shown in Figs.S13-S14.d,e, All forecast results for Spain and Germany.Shown are the RMSE forecast errors and their ratios.h, Average forecast errors for the 8 regions.The benefits of removing reporting delays vary by region.g, The average RMSE forecast error of experiment 2 as a function of SIRD model fitting error for all experiments.h, The mean relative forecast error as a function of SIRD model fitting error for all experiments.

TABLE I :
Naming convention for the epidemic data, Y = {I, R, D}.

TABLE II :
Inferred reporting delays.The table displays the inferred differences between the expected delay times E[T R ] − E[T D ], E[T I ] − E[T D ], standard deviations σ I , σ R and σ D , and optimized values for the objective function O b ( Ȳ ).See Tables