Understanding changes of laser backscattering imaging parameters through the kiwifruit softening process using time series analysis

ABSTRACT

During kiwifruit storage, quality monitoring is required for inventory planning and consistent quality maintenance. Commercial near-infrared (NIR) spectrophotometers showed reduced performance in the estimation of kiwifruit flesh firmness (FF) as the FF estimation is indirect and can be affected when both, textural structures and absorbing compounds, change during postharvest ripening. Laser backscattering imaging (LBI) records the backscattered signal after a single laser beam interacts with kiwifruit tissue, including merged information on light absorption and scattering. Measurements were carried out at 830 nm, where scattering is most dominant. In this work, time series of kiwifruit ‘Zesy002’ (n = 30) and ‘Hayward’ (n = 30) LBI were collected through the postharvest ripening during a 15-day shelf life at 20°C. Four LBI parameters capturing DIP, FWHM, SLP and HWQM were selected in this study. ‘‘Zesy002’ DIP, FWHM, SLP, and HWQM increased approx. 0.6 cm, 0.2 cm, 0.3 and 0.14 cm, respectively. ‘Hayward’ LBI increased approx. 0.2 cm, 0.1 cm, 0.2 and 0.04 cm, respectively. Different initial LBI values between cultivars and LBI changes may reflect the actual stage of softening, affected by kiwifruit ripeness. In conclusion, time series analysis could be useful in describing kiwifruit LBI change during ripening and making forecasts.

KEYWORDS:

• Actinidia
• Holt-winters
• ARIMA
• firmness
• forecast

Introduction

Kiwifruit ‘Zesy002’ (Actinidia chinensis var. chinensis) and ‘Hayward’ (A. chinensis var. deliciosa) are the two main commercial kiwifruit cultivars in New Zealand and worldwide. Kiwifruit is climacteric fruit, which is frequently harvested unripe and stored for long periods (6–8 months) in low temperature (0°C) storage conditions (Schröder and Atkinson 2006). This long-term storage allows kiwifruit to develop physiologically before becoming suitable for consumption (Beever and Hopkirk 1990). A three-phase softening curve has been reported for kiwifruit ripening, including a relatively stable initial lag phase, then a rapid exponential decline to 20–30 N, followed by a slow flesh firmness (FF) loss to eating ripeness (Hewett et al. 1999). The softening rate in the 2^nd phase and length of the 1^st phase may be influenced by multiple factors, including temperature, storage time, environmental ethylene concentration and maturity (Ritenour et al. 1999; Jabbar and East 2016; Burdon et al. 2017).

Monitoring the softening process during storage and throughout the logistics chain until fruit consumption is desirable for industry inventory decision making. The traditional destructive measurement of FF captures a small sub-sample from the entire batch of fruit, resulting in a batch estimation for the kiwifruit population. Several mechanistic models have been reported to quantify the kiwifruit FF softening trend. Particularly, a softening model was developed considering enzyme processes related to kiwifruit physiology of cell wall breakdown and ripening process (Hertog et al. 2016). The model had good performance with 97% of the data variation in FF explained. However, the feasibility of FF models is reduced, because the at-harvest FF has to be acquired destructively using compression firmness.

Non-destructive near-infrared (NIR) spectroscopy has been investigated for kiwifruit quality estimation, allowing the monitoring of quality changes of each individual or batch of fruit. NIR spectroscopy captures data containing the signal of absorption that is related to kiwifruit’s chemical properties. While some cell wall compounds may be detected, the signal also contains scattering effects. The challenge of FF measurement using NIR spectroscopy appears due to the varying stoichiometric relationship of FF and absorbing molecules during the dynamic softening process. Therefore, the indirect FF estimation based on NIR data is not consistently correlated to destructively measured FF levels (Walsh et al. 2020). Spatially-resolved spectroscopy (Lu et al. 2020) and time-resolved spectroscopy (Cubeddu et al. 2001) attempted to overcome the problem of working on a sum signal, by distinguishing the absorption and scattering part of the signal based on the discovered informative wavelengths from hyperspectral readings, however, relating optical measurements to the complex phenomenon of firmness remained challenging.

During kiwifruit softening, both cell wall composition and structural properties contribute to FF change. The rapid FF reduction results from pectin solubilisation where water-insoluble pectin becomes soluble (Redgwell and Percy 1992; Schröder and Atkinson 2006). Concurrent with pectin degradation, cells begin to swell and separate due to a loss of middle lamellae adhesion and eventually, the cell walls begin to break down (Redgwell et al. 1992). Besides cell wall degradation, other physical properties, such as water potential, osmotic potential, turgor pressure and tissue density, could also influence FF (Harker and Hallett 1994). During softening, kiwifruit turgidity decreases followed by an increase in air spaces and a decrease in water and osmotic potential (Harker and Hallett 1994). Ripe kiwifruit with the same turgor pressure as unripe kiwifruit were observed to have lower FF (Rojas et al. 2002).

Laser backscattering imaging (LBI) technology was developed to provide information on the photon attenuation profile while travelling through the fruit tissue (Lu et al. 2020). The radial attenuation profile could be analysed based on Farrell’s theory (Farrell et al. 1992) and optical properties, µ_a and µ_s’, are usually estimated via a curve-fitting process in the diffusion model. It has been frequently assumed that µ_s’ in the short-wave NIR range is correlated to kiwifruit FF, and µ_s’ could be estimated using LBI attenuation profile (Baranyai and Zude 2009). However, decoupling of µ_a and µ_s’ remains a major challenge for LBI when applied to fruit quality estimation (Baranyai and Zude 2009; Zude-Sasse et al. 2019). Alternatively, LBI parameters derived directly from the attenuation profile may be more feasible in identifying fruit and particularly kiwifruit with varying firmness as related to chilling injury symptoms (Sun et al. 2016; Yang et al. 2021).

LBI parameters may provide information on the changes in light attenuation during fruit ripening. The non-destructive character of the method allows repeated measurements of fruit during the ripening process. Time series analysis is a statistical method that builds a model to provide a plausible description of the observed data and facilitates trend forecasts on the sample development (Shumway and Stoffer 2017). The approach has been used in many fields, e.g. economy, and sport competition, but little work has been published on using time series to predict changes in fruit quality. Weng et al. (2019) applied an AutoRegressive Integrated Moving Average (ARIMA) model to predict agricultural product sales price and found it suitable to predict short periodical data rather than daily price estimation. Verroens et al. (2006) compared 4 time-series models to predict bell pepper cyclic production and found that utilising at harvest fruit size and hardness reached a better prediction of yield than using colour. González-Teruel et al. (2022) used a time series model to predict drought stress in cherry orchard based on soil moisture and provided an on-day alert for farmers. Times series analysis has the potential to investigate the correlation of data and decompose influences, thus interpreting the meaning of a data set.

This study aims to monitor the changes in LBI parameters that are associated with the softening of kiwifruit during shelf life ripening and attempts to use time series analysis methods for data interpretation. The specific objectives are to:

1. describe changes in kiwifruit LBI parameters during shelf life ripening;

2. determine whether LBI parameters are useful for describing kiwifruit softening patterns, and

3. identify possible kiwifruit quality changes associated with LBI parameters during ripening.

Materials and methods

Fruit source

‘Zesy002’ and ‘Hayward’ kiwifruit were commercially harvested from the Bay of Plenty region, New Zealand, in April and May 2021 respectively. Kiwifruit was graded and delivered to Massey University from commercial packhouses on 2^nd and 27^th May, respectively. ‘Zesy002’ was stored at 1°C for 26 days before being taken out of cool storage on 30^th May and warmed to room temperature (20°C) overnight prior to the start of the shelf-life study. ‘Hayward’ was stored at 0°C for 2 days with ripening starting simultaneously with ‘Zesy002’ after warming overnight on 31^st May. Fruit of uniform shape, and with no scar or sunburn observed on the surface were used.

Kiwifruit image and quality measurement

The LBI system (Figure 1A) included the light source with a laser diode (FP-D-DIG-520-17-C-F250-USB, Laser Components, Germany) emitting at 830 nm to minimise the influence on LBI parameters due to pigment absorption in kiwifruit flesh colour. The laser diode had 11.2 mW of power output. An 8-bit CCD camera (MER-131-210U3M NIR, China Daheng Group, China) providing a resolution of 1280 × 1024 pixels, equipped with F1.4 aperture zoom lens and 10–40 mm focal length (12VG1040 ASIR-SQ, Tamron Co. Ltd, Japan) was chosen. A desktop computer was used to control the light source and camera. In the geometric calibration, one pixel represented 0.1016$\times$0.1016 (mm²). The penetration depth can be increased by increasing laser power output or scattering (Jacques 2013). An electrical moving table was used to adjust the sample position to maintain a 25 cm distance between the highest point of the fruit surface and the camera. Kiwifruit LBI were taken in the equatorial region at 2 positions, 90° apart (Figure 1B).

Figure 1. A, Schematic of laser backscattering system; B, Kiwifruit LBI acquisition positions correspond to cross-marks.

The kiwifruit sample was segregated into 2 lots. Lot 1 kiwifruit was measured for the purpose of time series model development. Kiwifruit (n = 30) were selected as 3 fruits from each of 10 grower lines (GL) for each cultivar. The LBI parameters of each individual fruit were monitored for 13.5 days at 12-h intervals resulting in 27 LBI observations at 20°C. The initial acoustic firmness (AF) of each kiwifruit was measured using AWETA AFS (AWETA International Ltd., Nootdrop, Holland) at fruit shoulders around both stem and blossom ends, and the average reading of AFs from two locations was recorded as kiwifruit AF. The AF was selected to provide texture-related data on the same fruit measured with LBI.

Lot 2 kiwifruit (n = 40) from the same batch as lot 1 fruit was used as a reference for monitoring FF changes. For lot 2, 10 fruit per cultivar were measured on each of the 4 measurement occasions at 3-day intervals on days 1, 4, 7 and 10, representing the 1^st, 7^th, 13^th and 19^th LBI observations, resulting in a total of 40 fruit per cultivar. The FF was measured destructively after LBI acquisition by means of a compression test using a penetrometer (Willowbank Electronics Ltd., New Zealand) with a standard 7.9 mm diameter convex Effegi probe at a speed of 8 mm s⁻¹ to 8 mm depth.

Kiwifruit samples collected from lot 2 kiwifruit on days 1 and 10 were further subjected to microscopy following the same methodology as described previously (Li et al. 2021). Kiwifruit tissue samples were collected in the light diffusion area, corresponding to LBI image region between the solid orange circle and the solid blue circle in Figure 2, and stored in standard FFA fixative (100 ml fresh fixative contained 53 ml 95% ethanol, 10 ml 37% formaldehyde and 5 ml glacial acetic acid) at a 5°C chiller. A total of 20 kiwifruit samples per cultivar were fixed and processed. Kiwifruit samples were re-cut to fit and transferred to a wax embedding cassette. The cassettes with kiwifruit sections were dehydrated with a Leica TP1020 Semi-enclosed Benchtop Tissue Processor (Leica Microscopy Systems Ltd., Heerbrugg, Switzerland) with the process of dehydration in ethanol holding solutions at 50% (15 min), 70% (3 h), 95% (3 h), 100% (3 h), 100% (pure, 3 h) and 100% (pure, 3 h). For the wax infiltration, sections were processed with pure ethanol and xylene (1:1, 3 h), xylene (3 h) and xylene (3 h). The section was processed automatically over 2.5 days. Kiwifruit sections were then transferred from the processing basket to the wax hold bath. The embedding process used 3.0 cm $\times$ 2.4 cm mould with a kiwifruit sample sitting in the middle using a Leica EG1160 paraffin dispensing system (Leica Biosystems, Heerbrugg, Switzerland). Kiwifruit sections were orientated with the outer pericarp side facing the short sides of the embedding cassettes. Once the samples were stable, they were stored at room temperature. Embedded samples were sectioned using a Leitz 1512 rotary microtome (Leitz, Germany) and a 0.25 mm disposable microtome steel blade. The wax samples were trimmed into a triangular shape first. The 15 µm samples were sectioned, placed in a 42°C water bath (Leica HI1210, Leica Biosystems, Heerbrugg, Switzerland) then collected onto a slide and dried at 40°C overnight. Dry sections on slides were stained in the Coplin staining jar with a 0.05% aqueous solution of toluidine blue for 5 minutes. Sections then were rinsed in water and left to dry for 30 minutes at 40°C then deparaffinised. Sections were taken through xylene 1–4 for 15 minutes each in Coplin staining jars, then air-dried and permanently mounted with a 24 $\times$ 50 mm glass coverslip using TBS SHUR/Mount™ liquid mounting medium (Triangle Biomedical Science, North Carolina, United States). Mounted slides were left in a fume hood overnight at room temperature before microscopic evaluation. The samples were analysed using an Olympus light microscope (Vanox AHT3, Olympus Optical Co. Ltd., Tokyo, Japan) and the images were captured and visualised by a colour camera (DP74, Olympus Optical Co. Ltd., Tokyo, Japan) and cellSens standard digital imaging software version 1.7 (Olympus Optical Co. Ltd., Tokyo, Japan).

Figure 2. Laser backscattering image (LBI) and attenuation profile. A, Raw image obtained from a kiwifruit sample at 830 nm with an image size of 65.02 $\times$ 52.02 (mm²). B, Transformed diffusion profile (attenuation profile) from LBI. LBI parameters were the radius of the saturated area or distance to incident point (DIP), the radius at 75% of maximum intensity (Q1R), full width at 50% of maximum intensity (FWHM), the radius at 25% of maximum intensity (Q3R), and the slope of the linear regression model (SLP) built with log-transformed profile data between Q1R and Q3R. The blue cross in A represents the incident point, the orange circle in solid line in A is the saturated area with the radius corresponding to DIP in B, the blue circle in A is around at 1 cm radius B, and apparent light diffusion stops at the orange dashed circle in A.

Time series analysis input parameters

LBI images were analysed in R (version 3.6.0, R Foundation for Statistical Computing, Vienna) running in RStudio (version 1.0.153, RStudio Inc.). The raw image was loaded and normalised into a grey-scale matrix using the ‘readbitmap’ package (version 0.1.5). After the normalisation process, light intensity at each pixel has a value between 0 and 1. Normalisation could determine the saturated area first and then calculate the weighted centre. The resulting grey-scale matrix contains normalised intensity and coordinates. The saturated area has a light intensity of 1, and the light intensity reduces to 0 at the end of light diffusion (Figure 2A). The coordinates of the weighted centre of the saturated area are calculated as incident point coordinates transformed by geometric calibration in cm. Pixels at the same distance to the incident point can be calculated in 1-pixel concentric rings around the incident point. The average backscattering intensity of the ring at each radius was computed into LBI light attenuation profiles (Figure 2B).

Extracted parameters (Figure 2B) of the profile were full width at half maximum intensity (FWHM, Eq. 2.1), the radial distance of 75% (Q1R, Eq. 2.2) and 25% intensity values (Q3R, Eq. 2.3), the half width of quarter maximum intensity (HWQM) was calculated using FWHM and Q3R to describe the radius between 25% – 50% of the maximum intensity (Eq. 2.4). and the slope of the linear regression model (SLP) built with log-transformed profile data between Q1R and Q3R. 2.1 $\begin{array}{r}F\end{array}\mathrm{WHM}=2\cdot \left(r_{I_{m50}}-r_{I_{\max }}\right)$2.1

2.2 $Q1R=r_{I_{m75}}-r_{I_{\max }}$2.2 2.3 $Q3R=r_{I_{m25}}-r_{I_{\max }}$2.3 2.4 $\mathrm{HWQM}=Q3R\text{–}\left({\displaystyle \frac{\mathrm{FWHM}}{2}}\right)$2.4 With r = radial distance measured from the incident point, $I_{\max }$ = maximum intensity value (corresponding to the orange solid circle area in Figure 2A, and plateau phase in Figure 2B), $I_{m75}$ = 75% of maximum intensity value (Eq. 2.5), $I_{m50}$ = 50% of maximum intensity value (Eq. 2.6), $I_{m25}$ = 25% of maximum intensity value $I_{\min }$ = minimum intensity value (Eq. 2.7). 2.5 $I_{m75}=I_{\max }-0.25\left(I_{\max }-I_{\min }\right)$2.5 2.6 $I_{m50}=I_{\max }-0.5\left(I_{\max }-I_{\min }\right)$2.6 2.7 $I_{m25}=I_{\min }+0.25\left(I_{\max }-I_{\min }\right)$2.7 The average LBI parameters from two positions of each kiwifruit were used for time series analysis model development and forecasting. The time series model of each LBI parameter was developed using a simple linear regression model, Holt-Winters smoothing (HWS) and AutoRegressive Integrated Moving Average (ARIMA). Details are provided in the following sections.

Time series analysis models of LBI parameters were developed based on the first 20 observations, and the model performance assessment was conducted on the last 7 values forecasted using the developed models and compared with the real observations. Fitted values using a regression model were calculated with Eq. 2.8. The forecasting via HWS and ARIMA was conducted using the ‘Forecast’ package (version 8.15). Data smoothing methods were applied to the averaged LBI data in HWS and ARIMA for comparing the performance in the forecasting.

Time series analysis with linear regression model

The regression of LBI parameter and time model was built by (Eq. 2.8).

2.8 $\widehat{y_t}=a{t}^2+\mathrm{bt}+c$2.8

Where $\widehat{y_t}$ is the fitted value at t^th observation, each observation was at 12-h time interval. a and b are polynomial coefficients, c is a constant. The model is fitted with the average LBI parameters of 30 fruit.

In the regression model evaluation, the mean of residuals $e_t$ (Eq. 2.9) is targeting zero with the lowest minimum sum of squared errors (SSE, Eq. 2.10) and maximum adjusted R² (${\overline{R}}^2$, Eq. 2.11). 2.9 $e_t=y_t-\widehat{y_t}$2.9 2.10 $\mathrm{SSE}=\sum _{t=1}^N{e}_t^2$2.10 2.11 ${\overline{R}}^2=1-\left(1-R^2\right){\displaystyle \frac{N-1}{N-k-1}}$2.11

Where $y_t$ is the observed value and $\widehat{y_t}$ is the fitted value at t^th observation. $e_t$ is the residual at t^th observation. N is the total number of observations t and k is the number of predictors (k = 3 in this study).

Time series analysis with Holt-Winters smoothing

Holt-Winters smoothing (HWS) forecasting is described below in Eq. 2.12 – Eq. 2.15 (Holt 2004): 2.12 $S_t=\alpha \cdot y_t/\left(I_t\text{–}L\right)+\left(1\text{–}\alpha \right)\left(S_{t-1}+b_{t-1}\right)$2.12 2.13 $b_t=\text{β}\cdot \left(S_t\text{–}{S}_{t-1}\right)+\left(1\text{–}\text{β}\right)\cdot b_{t-1}$2.13 2.14 $I_t=\gamma \cdot y_t/S_t+\left(1\text{–}\gamma \right)I_{t-\mathrm{L}}$2.14 2.15 $F_{t+m}=\left(S_t+m{b}_{t-1}\right)\cdot I_{t-L+m}$2.15 Where α is the level smoothing factor, β is the trend smoothing factor, γ is the seasonality smoothing factor, y_t is the actual value at t^th observation, S_t is the smoothed value at t^th observation, b_t is the best linear trend regardless of seasonal changes and I_t is the seasonal factor, F_t + m is the forecast value at m steps ahead from t^th observation, L is cycle length (Holt 2004).

Unlike regression with smoothing which allocates the same weight to both the previous and the newly observed value, the HWS exponentially varies weights and values towards recent data (Gelper et al. 2010). In HWS, the two parameters α and β represent the weight of the recent LBI parameter data and its trend, respectively. The range of α and β are between 0 and 1. In this study, no seasonality or cyclic pattern was evaluated and, therefore, no third coefficient was introduced. In kiwifruit postharvest softening, a single decreasing trend is expected. Holt (2004) explained that a smaller α and β (closer to 0) indicate a greater smoothing effect, which means that forecasting is based on a longer period of the previous level (L) and a long-term trend (T), respectively. Alternatively, a larger α and β (closer to 1) indicates less smoothing effects with the forecasting giving greater weight to the recent L and T in the data. When α and β equal 1, the forecast output is the same as the last actual observation.

Time series analysis with ARIMA

ARIMA (q, d, p) is usually employed for investigating the irregular component of a time series (Box et al. 2015). The general forecasting process introduced by Box et al. (2015) using ARIMA includes firstly finding the differencing coefficient (d, Eq. 2.16 – Eq. 2.17) to make the data set stationary. Once the distribution of observation ($y_t$) is not affected by the time points, the data is stationary (Kwiatkowski et al. 1992). Secondly, a residual test is conducted to determine the order for autoregressive (AR, p, Eq. 2.18) and moving average (MA, q, Eq. 2.19) by autocorrelation function (ACF) and partial autocorrelation function (PACF) tests. Thirdly, the models are compared using Akaike information criterion (AIC, Eq. 2.30) and Bayesian information criterion (BIC, Eq. 2.31) to identify the best model with the minimum AIC and BIC.

The d or integration (I) in ARIMA model computes the differences between observations to make non-stationary time series stationary to exclude the underlying correlation with time. The first-order differenced time series (d = 1) can be written as Eq. 2.16. 2.16 $y{{}^{\mathrm{\prime }}}_t=y_t-y_{t-1}$2.16 Where $y_t$ is the observation at t^th observation point with d = 1 and t ≥ 2. $y_t^{\mathrm{\prime }}$ is the observation at t^th after the differentiation.

The second-order differenced time series (d = 2) can be written as Eq. 2.17. 2.17 $y{}^{\mathrm{\prime }}{{}^{\mathrm{\prime }}}_t=y_t-2y_{t-1}-y_{t-2}$2.17 Where $y_t$ is the observation at t^th observation point with d = 2 and t ≥ 3. $y{}^{\mathrm{\prime }}{{}^{\mathrm{\prime }}}_t$ is the observation at t^th after the second differentiation.

The order of d refers to the successive first differences, and $y_t^{\mathrm{\prime }}$ represents the mean of time series by removing the change of level thus the trend influence can be limited (Kwiatkowski et al. 1992). In general, when there is constant c (Eq. 2.18 – Eq. 2.20), the forecasted values could be the mean of data (when d = 0), following a straight line (when d = 1) or in a quadratic trend (when d = 2). In this study, kiwifruit FF is expected to have d = 2 as an exponentially decreasing FF trend was expected.

The AR (p) in ARIMA model utilises a linear combination of the past values to forecast future values as described in Eq. 2.18. The MA (q) in ARIMA model uses past forecast errors to make forecasts (Eq. 2.19). In this study of kiwifruit ripening, non-seasonal ARIMA is expected and can be written as Eq. 2.20. 2.18 $y_t=c+{\mathrm{\varnothing }}_1{y}_{t-1}+{\mathrm{\varnothing }}_2{y}_{t-2}+\cdots +{\mathrm{\varnothing }}_p{y}_{t-p}+{ϵ}_t$2.18 2.19 $y_t=c+{ϵ}_t+{\theta }_1{ϵ}_{t-1}+{\theta }_2{ϵ}_{t-1}+\cdots +{\theta }_q{ϵ}_{t-q}$2.19 2.20 $y{{}^{\mathrm{\prime }}}_t=c+{\mathrm{\varnothing }}_{1}y{{}^{\mathrm{\prime }}}_{t-1}+{\mathrm{\varnothing }}_{2}y{{}^{\mathrm{\prime }}}_{t-2}+\cdots +{\mathrm{\varnothing }}_{p}y{{}^{\mathrm{\prime }}}_{t-p}+{\theta }_1{ϵ}_{t-1}+{\theta }_2{ϵ}_{t-1}+\cdots +{\theta }_q{ϵ}_{t-q}+{ϵ}_t$2.20 Where ${ϵ}_t$ is the white noise, ${\mathrm{\varnothing }}_1$, … , ${\mathrm{\varnothing }}_p$ are parameters of the lagged values of $y_t$. These parameters determine the pattern of the time series data. The restrictions of the value of $\mathrm{\varnothing }$ for AR (p = 1) are Eq. 2.21 and for AR (p = 2) are Eq. 2.22-Eq. 2.24. Similarly, the restrictions of the value of $\theta$ for MA (q = 1) are Eq. 2.25 and for MA (q = 2) are Eq. 2.26-Eq. 2.28 (Box et al. 2015). 2.21 $-1<{\mathrm{\varnothing }}_1<1$2.21 2.22 $-1<{\mathrm{\varnothing }}_2<1$2.22 2.23 ${\mathrm{\varnothing }}_1+{\mathrm{\varnothing }}_2<1$2.23 2.24 ${\mathrm{\varnothing }}_2-{\mathrm{\varnothing }}_1<1$2.24 2.25 $-1<{\theta }_1<1$2.25 2.26 $-1<{\theta }_2<1$2.26 2.27 ${\theta }_1+{\theta }_2<1$2.27 2.28 ${\theta }_1-{\theta }_2<1$2.28

The value of AR order (p) and the MA order (q) could be considered using ACF and PACF plots comparing the pattern (either gradually decreasing or significant cut at certain lags). The maximum likelihood estimation (MLE), AIC and BIC are used for model evaluation. The best model will be selected with the maximum MLE and minimum AIC and BIC (Eq. 2.29 – Eq. 2.31). 2.29 $\mathrm{MLE}=\sum _{t=1}^N{ϵ}_t^2$2.29 2.30 $\mathrm{AIC}=\mathrm{Nlog}\left({\displaystyle \frac{\mathrm{SSE}}{N}}\right)+2\left(k+2\right)$2.30 2.31 $\mathrm{BIC}=\mathrm{Nlog}\left({\displaystyle \frac{\mathrm{SSE}}{N}}\right)+\left(k+2\right)\log \left(N\right)$2.31 where N is the total number of observations used for estimation and k is the number of predictors in the model.

Results and discussion

LBI parameter analysis

LBI parameter trends

The DIP represents the saturated area around the incident point of the laser beam in light diffusion inside kiwifruit. A higher DIP indicates that the photon detectable signal started to reduce at a larger radius from the incident point (Figure 2). The SLP represents the rate of attenuation of backscattering intensity between 75% and 25% of maximum intensity in the diffusion region. Hence higher SLP indicates a longer travel path of photons. The ½ FWHM and HWQM represent the diffusion area within 75% to 50% and 50% to 25%, respectively. Higher FWHM and HWQM also indicate less photon energy reduction at the first and the second half photon travel region respectively.

‘Zesy002’ kiwifruit DIP, SLP and HWQM showed a quadratic increasing trend over the 27 observations captured (Figure 3). The increased LBI parameters of DIP (Figure 3A), SLP (Figure 3C) and HWQM (Figure 3D) indicate that the effective photon’s travel path became longer during storage at 20°C. The average LBI parameter of lot 2 kiwifruit (Figure 3 orange dashed lines) showed similar trends as lot 1 kiwifruit (Figure 3 black solid lines), and FF decreased during storage (Figure 3 blue dotdashed lines). Therefore, FF reduction may increase LBI parameter values resulting from physical property changes. However, a large variation of LBI parameters was observed for different fruit samples on each measurement day in lot 2. Coincidingly, FF variation was observed on day 1 and day 4. Kiwifruit used for time series analysis in lot 1 had less variation of LBI parameters compared with lot 2 fruit samples, as the same kiwifruit samples were measured on each measurement occasion in that case. The increase of FWHM in the later storage may indicate that FWHM may be sensitive at the lower FF range and could determine the stage of the softening process when FF reduction plateaus. In the previous work (Yang et al. 2021), low correlation was observed between Q1R, Q3R and kiwifruit with chilling injury symptoms, thus those two LBI parameters were excluded in the current work.

Figure 3. Changes of laser backscattering imaging (LBI) profile parameters at 830 nm during 14 days at 20°C for ‘Zesy002’ kiwifruit. LBI parameters are DIP (A, the radius of the saturated area), FWHM (B, double the radius at 50% of maximum intensity), SLP (c, the slope of the linear regression model built with log-transformed profile data between 25% and 75% of maximum intensity) and HWQM (d, the radius between 25% and 50% of maximum intensity). Coloured solid lines are the average LBI parameter for each grower line (n = 3), the black solid line is the average LBI parameter for all grower lines (n = 30). The blue dotdashed line is the average firmness and the orange dashed line is the average LBI parameter for kiwifruit (n = 10) sourced from the same batch. Error bar is standard deviation.

A similar result was observed for ‘Hayward’ kiwifruit (Figure 4), all LBI parameters showed an increasing trend through ripening. In addition, a stable lag phase of ‘Hayward’ was observed from day 1 to day 6 for DIP, SLP and HWQM (Figure 4A,C,D), whilst a longer lag phase was observed for FWHM (Figure 4B). No lag phase after the increasing period was observed for all LBI parameters. The negative correlation of LBI parameters and FF was also observed for ‘Hayward’, where FF variation was reduced on day 10.

Simple regression model

For the first 20 observations of ‘Zesy002’ kiwifruit, the regression model fitted well for DIP, SLP and HWQM with a high R² over 0.99 (Table 1). The model has a positive coefficient a, reflecting the similar pattern between LBI parameters and storage time. A lower R² was obtained for FWHM where no clear trend was observed (Figure 3B). An interpretation of c is the average initial value of LBI parameters when t = 0.

Figure 4. Changes of laser backscattering imaging (LBI) profile parameters at 830 nm during 14 days at 20°C for ‘Hayward’ kiwifruit. LBI parameters are DIP (A, the radius of the saturated area), FWHM (B, double the radius at 50% of maximum intensity), SLP (c, the slope of the linear regression model built with log-transformed profile data between 25% and 75% of maximum intensity) and HWQM (d, the radius between 25% and 50% of maximum intensity). Coloured solid lines are average LBI parameter for each grower line (n = 3), the black solid line is average LBI parameter for all grower lines (n = 30). The blue dotdashed line is average firmness and orange dashed line is average LBI parameter for kiwifruit (n = 10) sourced from the same batch. Error bar is standard deviation.

Table 1. Kiwifruit laser backscattering imaging (LBI) parameter analysis with fitted valued at each observation data point using simple linear regression.

LBI	a	b	c	R^2	SSE
‘Zesy002’
DIP	1.05e^−03	−8.27e^−05	2.17	0.99	2.76e^−03
FWHM	5.33e^−04	−6.53e^−03	2.03	0.77	5.79e^−03
SLP	3.35e^−04	4.08e^−03	−1.12	0.99	8.87e^−04
HWQM	1.24e^−04	8.11e^−04	2.67	0.99	5.10e^−05
‘Hayward’
DIP	4.29e^−04	−1.75e^−03	1.96	0.97	1.44e^−03
FWHM	8.98e^−05	−7.46e^−05	1.91E	0.64	1.31e^−03
SLP	3.52e^−04	−1.79e^−03	−1.32	0.97	6.90e^−04
HWQM	5.65e^−05	−1.05e^−04	2.22	0.95	4.69e^−05

For ‘Hayward’ kiwifruit, the R² of the simple regression model was lower than for ‘Zesy002’ with the R² being over 0.9 in DIP, SLP and HWQM (Table 1). A low R² was also observed for ‘Hayward’ FWHM where no clear trend was observed (Figure 4B). The coefficients of the simple regression model describe the general trends based on previous observation data, however, without identifying the phase switching point during the kiwifruit softening/ripening process. Putting the same weight on each observation may not well represent future changes during ripening.

Holt-Winters smoothing model

In the HWS model (Table 2), DIP and FWHM for ‘Zesy002’ had a lower α (0.3 and 0.1) but higher β (1.0), which indicates that the forecasted L is based on long-term observation, while the forecasted T is built on the most recent data. HWQM showed an opposite pattern (α = 0.7, β = 0.3), indicating that the forecasted L was based on recent data while the T was on a long-term basis. For SLP, both α and β were approx 0.5, thus the L and T forecasts were equally influenced by the past data regardless of the observation period. ‘Zesy002’ SLP did steadily increase with no lag phase observed (Figure 3C), thus the forecasts didn’t put particular weight on historical data or more recent data. On the other hand, a lag phase in DIP (Figure 3A) and FWHM (Figure 3B) was observed followed by a period of increase, where a stable forecast is expected (α = 0.3 and 0.1). For ‘Zesy002’ HWQM (Figure 3D), no lag phase was observed after a steady increase. Therefore, the forecasted value will follow the recent T, or increasing rate, and continue to increase (α = 0.7).

Table 2. Kiwifruit laser backscattering imaging (LBI) parameter analysis with fitted valued at each observation data point using Holt-Winters Smoothing.

LBI	a	b	α	β	SSE
‘Zesy002’
DIP	2.58	0.03	0.33	1.00	0.01
FWHM	2.11	0.02	0.12	1.00	0.01
SLP	−0.91	0.01	0.51	0.42	0.00
HWQM	0.33	0.01	0.73	0.30	0.00
‘Hayward’
DIP	2.09	0.01	0.83	0.00	0.00
FWHM	1.94	0.00	0.61	0.29	0.00
SLP	−1.21	0.01	0.82	0.78	0.00
HWQM	0.24	0.00	1.00	0.28	0.00

In ‘Hayward’ HWS model (Table 2), DIP and FWHM had a short-term L forecast (α = 0.8 and 0.6) and a long-term T forecast (β = 0 and 0.3). This finding for ‘Hayward’ DIP and FWHM was opposite with ‘Zesy002’. DIP (Figure 4A) and FWHM (Figure 4B) had a lag phase prior to an increasing phase. HWQM (Table 2) showed the same short-term L forecast and long-term T forecast (α = 1.0, β = 0.3). SLP (Table 2) has both high α (0.8) and β (0.8), indicating that L and T will follow the recent data, and are expected to continue increasing.

ARIMA model

The selected ARIMA models of ‘Zesy002’ LBI parameter are presented. DIP, SLP and HWQM had a d = 2, which indicated a quadratic trend of forecasting. Whilst FWHM had a d = 1, indicating a straight line of forecasted values. The d values of each parameter agree with the observation in Figure 3. HWQM had the highest log-likelihood and lowest AIC, hence better performance of forecasting is expected. DIP and FWHM had a lower log-likelihood and higher AIC, thus the ARIMA model might not be suitable for analysing ‘Zesy002’ ripening data (Table 3).

Table 3. Kiwifruit laser backscattering imaging (LBI) parameter analysis with fitted value at each observation data point using Auto Regressive Integrated Moving Average (ARIMA) model.

LBI	p	d	q	AIC	loglik
‘Zesy002’
DIP	1.00	2.00	2.00	−87.87	47.93
FWHM	1.00	1.00	0.00	−84.29	44.15
SLP	0.00	2.00	1.00	−110.65	57.33
HWQM	0.00	2.00	1.00	−156.36	80.18
‘Hayward’
DIP	0.00	1.00	0.00	−108.31	55.15
FWHM	1.00	1.00	0.00	−115.48	59.74
SLP	1.00	2.00	0.00	−131.67	67.83
HWQM	0.00	2.00	0.00	−175.80	88.90

The selected ARIMA model of LBI parameters had a higher log-likelihood and lower AIC for ‘Hayward’ HWQM and SLP. Alternatively, less variation was observed in SLP and HWQM than DIP and FWHM (Figure 4). Looking at d values, SLP and HWQM were described similar to ‘Zesy002’ with a quadric trend. The straight trend of DIP (d = 1) could be due to missing the second lag phase (Figure 4A).

For all LBI parameters, residuals of fitted values were within the 95% confidence interval (Figure 5). Residuals of the later stage forecasted values were outside of the threshold of 95% confidence interval.

Figure 5. Residuals analysis of times series analysis models. Residuals of fitted (empty) and forecasted (closed) values using Simple regression, Holt-Winters smoothing (HWS) and Auto Regressive Integrated Moving Average (ARIMA) model were compared with measured data. LBI parameters of DIP (A, the radius of the saturated area), FWHM (B, double the radius at 50% of maximum intensity), SLP (C, the slope of the linear regression model built with log-transformed profile data between 25% and 75% of maximum intensity) and HWQM (D, the radius between 25% and 50% of maximum intensity) are presented. Dash lines are the upper and lower residual thresholds for a 95% confidence interval.

LBI parameter forecasting using time series analysis

For both kiwifruit cultivars, the forecast of the final 7 observations showed that regression was the best model to forecast HWQM, while HWS and ARIMA provided the lowest deviation for the prediction of DIP and FWHM (Figure 6). Although neither HWS nor ARIMA achieved good fitted results (first 20 observations) and forecasted results (last 7 observations) for FWHM (Figure 5), HWS still obtained some fluctuation during the ripening process. The ARIMA model is not ideal for kiwifruit shelf-life ripening data analysis due to the data size, as the data set is reduced after differencing and smoothing (MA process), which could potentially loose key features in the model. On the other hand, ARIMA is more powerful when more data becoming available in automated LBI measurements.

Figure 6. Fitted and forecasted kiwifruit laser backscattering imaging (LBI) profile parameters. DIP (A,B, the radius of the saturated area), FWHM (C,D, double the radius at 50% of maximum intensity), SLP (E,F, the slope of the linear regression model built with log-transformed profile data between 25% and 75% of maximum intensity) and HWQM (G,H, the radius between 25% and 50% of maximum intensity) for “Zesy002’ (A,C,E,G) and ‘Hayward’ kiwifruit (B,D,F,H) for 14 days at 20°C measured at 830 nm. Simple regression, Holt-Winters smoothing (HWS) and Auto Regressive Integrated Moving Average (ARIMA) model were compared with measured data.

Forecasted DIP was overestimated for ‘Zesy002’ by 3% (Regression) and 1% (HWS and ARIMA)and underestimated for ‘Hayward’ by 1% (HWS) and 3% (ARIMA, Figure 5A, Figure 6A,B). For ‘Zesy002’, the forecasted FWHM were underestimated using ARIMA model (2%) but were overestimated using regression or HWS in contrast (1%, Figure 5B, Figure 6C). For ‘Hayward’ the forecasted FWHM was underestimated using all time series models (1%, Figure 5B, Figure 6D). All forecasted SLP residuals were within the 95% confidence interval for both cultivars except regression for ‘Hayward’ (Figure 5C), with most SLP being underestimated for ‘Zesy002’ by 2% (ARIMA, HWS) and overestimated for ‘Hayward’ by 1% (ARIMA, HWS, Figure 6E,F). Forecasted ‘Hayward’ HWQM residuals were out of 95% confidence interval using ARIMA (Figure 5D), and forecasted HWQM were underestimated for both cultivars by 2% (Figure 6G,H).

Discussion of LBI parameter changes associated with kiwifruit ripening

LBI parameter and firmness

‘Zesy002’ log-transformed HWQM showed the highest correlation to log-transformed FF (0.7), and ‘Hayward’ log-transformed negative SLP showed the highest correlation to log-transformed FF (0.8). For both cultivars, parameter FWHM had the lowest correlation (Table 4). In general, LBI parameters had a higher correlation to FF for ‘Hayward’ than ‘Zesy002’, which may be due to ‘Hayward’ having a wider FF range. During the experimental period, FF reached 3.7 N (from 57 N) and 2.3 N (from 32 N) for ‘Hayward’ and ‘Zesy002’, respectively.

Table 4. Correlation of, both log-transformed, LBI parameters and kiwifruit flesh firmness for ‘Zesy002’ and ‘Hayward’.

	DIP	FWHM	SLP	HWQM
‘Zesy002’	−0.65	−0.11	−0.63	−0.68
‘Hayward’	−0.75	−0.37	−0.78	−0.70

Numerical values are Pearson correlation coefficients. LBI parameters at 830 nm include the radius of the saturated area (DIP), double the radius at 50% of maximum intensity (FWHM), and the slope of the linear regression model (SLP) built with log-transformed profile data between the radius at 75% of maximum intensity and the radius at 25% of maximum intensity, the radius between 50% and 25% of maximum intensity (HWQM).

For both cultivars, higher LBI parameters were observed with lower FF of lot 2 kiwifruit in the early storage stages (Figure 7). A large variation of LBI parameters was observed when FF was over 15 N, however, a separation tendency of firmer fruit with lower LBI parameters and softer fruit with higher LBI parameters was observed (Figure 7).

Figure 7. Laser backscattering imaging (LBI) profile parameters and kiwifruit firmness correlation over 10 days at 20°C for ‘Zesy002’ kiwifruit (red) and ‘Hayward’ kiwifruit (green) at 830 nm. LBI parameters are DIP (A, the radius of the saturated area), FWHM (B, double the radius at 50% of maximum intensity), SLP (C, the slope of the linear regression model built with log-transformed profile data between 25% and 75% of maximum intensity) and HWQM (D, the radius between 25% and 50% of maximum intensity).

Combined LBI parameters (DIP, Q1R, FWHM, Q3R, SLP) demonstrated the ability to segregate unmarketable kiwifruit with chilling injury and resulted in 92% true positive cases in ‘Zesy002’ and 39% true positive cases in ‘Hayward’ (Yang et al. 2021). In future studies, the forecasted LBI parameters in ‘Zesy002’ can be tested in discriminant analysis of firmness classes by means of a similar approach. However, a different experiment requiring fewer points in time, but more fruit in shelf life for validation will need to be carried out.

LBI parameter and microstructure

The increase in LBI parameters along with the decrease in FF (Figures 3, 4, 7) could be explained by how microstructure influences optical properties. During kiwifruit softening, kiwifruit microstructural changes include cell wall swelling and pectin disintegration (Hallett et al. 1992). Cell swelling, cell degradation and shrivelling could be observed (Figure 8B,D) on day 10 for both cultivars, which could affect LBI parameters. For ‘Hayward’ kiwifruit, bruised fruit had higher μ_a and lower μ_s’ during storage ripening due to the change of tissue colour, FF and microstructure (Gao et al. 2021). Besides, microstructure information obtained from X-ray CT was used for studying fruit ripening processes and linked with scattering properties in mango (Cantre et al. 2014) and apple (Wang et al., 2020b).

Figure 8. Kiwifruit microstructure image on day 1 (A,C) and day 10 (B,D) at 20°C for ‘Zesy002’ kiwifruit (A,B) and ‘Hayward’ (C,D). Kiwifruit microstructure in the subsurface region and bar = 500 µm.

Different LBI parameters are associated with optical properties. The different regions of the LBI profile (Figure 2), defined as LBI parameters, possibly carry information on absorption and scattering in spatially-resolved system (Cubeddu et al. 1996). Limit work was reported on kiwifruit LBI parameters, while research has been conducted on kiwifruit optical properties. Low FF resulting from cell wall degradation followed the same pattern of increasing LBI parameters (Figures 3 and 4), while high LBI parameters could be a result of either low µ_a or low µ_s’ or both (Baranyai and Zude 2009). Destructive inverse adding double method, ‘Zesy002’ and ‘Hayward’ data of µ_s’ (cm⁻¹) ranged between 1 and 3 and 3 and 10 at 400–100 nm and 900–1700nm, respectively (Gao et al. 2021; Wang et al., 2020a). Another green kiwifruit ‘Huayou’ had higher µ_s’ (6 and 16 cm⁻¹) at the same NIR band (Liu et al. 2019). Non-destructive µ_s’ results were 13 and 16 cm⁻¹ using a time-resolved system.

Discussion of kiwifruit LBI parameter trends using Holt-Winter smoothing analysis

Absolute LBI parameter values of ‘Hayward’ kiwifruit were lower (p < 0.01) compared with ‘Zesy002’ throughout the experimental period (Figure 6). Besides optical properties that may differ between kiwifruit cultivars, the lower LBI may be due to ‘Zesy002’ and ‘Hayward’ measured in different ripening stages. ‘Zesy002’ fruit was in an advanced stage at the beginning of time series measurement. This assumption is consistent with the AF measurement of that ‘Zesy002’ AF (18 10⁶ Hz² g^2/3) being lower at the beginning of the experiment than ‘Hayward’ AF (27 10⁶ Hz² g^2/3), which may be a result of ‘Zesy002’ being stored for 26 days prior to measurement. The initial AF of ‘Zesy002’ and ‘Hayward’ is equivalent to approx. 32 and 57 N FF (Yang, unpublished data). Non-destructive AF was reported with a high segmented regression R² (0.89) within the same season, however, a large variation was observed when kiwifruit FF ranged between 20 and 60 N (Sneddon et al. 2024). Therefore, AWETA AF wouldn’t be suitable for the FF range interested in this study and for the sorting line adaption. Coordinately, the initial LBI parameters of ‘Zesy002’ (2.18 cm, – 1.1, and 0.27 cm for DIP, SLP, and HWQM, respectively) were higher than ‘Hayward’ (2.10 cm, – 1.3, and 0.25 cm for DIP, SLP, and HWQM, respectively). During the experimental period, ‘Hayward’ also maintained consistently higher FF (3.7 N at day 20) compared with ‘Zesy002’ (2.3 N at day 20). ‘Hayward’ LBI parameters on day 13 (2.20 cm, 1.15 and 0.26 cm for DIP, SLP and HWQM, respectively) were similar to ‘Zesy002’ on day 1.

The trend of LBI parameters may vary at different ripening stages, which corresponds to the β value in HWS modelling (Table 2). The fruit condition of ‘Zesy002’ was assumed within the 2^nd and 3^rd phases of the softening process, which captured a rapid decreasing and slowly decreasing trend, resulting in a lower FF decreasing rate in the 3^rd phase. Coordinately, SLP and HWQM, which have a high correlation with FF (Table 4), had a lower β in HWS analysis for ‘Zesy002’. The lower β indicated that the increasing rate of SLP and HWQM may become lower, when SLP and HWQM may pass the rapidly increasing period and entered the plateau stage. This is also supported by residuals of forecasted SLP and HWQM being within a 95% confidence interval (Figure 5C,D). The fruit condition of ‘Hayward’ was assumed within the 1^st and 2^nd phases of the softening process (Burdon et al. 2017), which captured a stable lag phase and a rapid decrease phase of FF. Coordinately, a stable lag phase was observed from day 1 to day 6 for ‘Hayward (Figure 4), whilst no lag phase was observed for ‘Zesy002’ (Figure 3). Since the LBI data for time series analysis may be collected when FF was within 1^st and 2^nd softening phases, whilst forecasts may be carried out about to enter or just entered the 3^rd softening phase i.e. slower decreasing FF phase, the forecasted LBI parameters were likely to be over-estimated (Table 2). The forecasted SLP residuals for ‘Hayward’ supported this assumption (Figure 5C). ‘Hayward’ DIP and HWQM in HWS analysis had a higher α but lower β (Table 2), which indicated that LBI parameters were continuously increasing. However, the increasing trend effect was diluted due to the lag phase, thus resulting in an underestimation of forecasted DIP and HWQM (Figure 5A,D).

The α and β were similar (with lower β than α) for SLP in HWS analysis for both ‘Zesy002’ and ‘Hayward’. Therefore, the value of α and β of LBI might be a single useful information that represents the ripeness. Lower α and β of ‘Zesy002’ SLP may be due to advanced ripening than ‘Hayward’. For ‘Hayward’ SLP, higher α and β determined the forecast L and T. Therefore, instead of forecasting LBI parameters, the α and β of SLP may be useful to analyse which phase the kiwifruit is located in. Apart from α and β in HWS, there might be other potential LBI parameters. For both cultivars, DIP, FWHM and HWQM followed the same trend but showed varied phase lengths (for example. HWQM entered the 2^nd lag phase the last), the ratio of LBI parameters might be useful in determining the ripening stage instead of single LBI parameters.

Although the starting point of the rapidly increasing phase in LBI parameter curve could be influenced by the at-harvest maturity of kiwifruit. In this study, both cultivars were selected from the main harvest. Thus each cultivar was assumed to have similar at harvest maturity and the miss of the initial lag phase of ‘Zesy002’ was assumed due to the extra postharvest storage time.

Discussion of time series model performance

Time series model performance could be influenced by an incomplete data set used for fitting. Both ‘Zesy002’ and ‘Hayward’ LBI parameters curves may be incomplete due to kiwifruit softening in this work were not through all 3 phases. Therefore, Figure 3 might represent the ‘Middle to end’ (5–60 N of FF) LBI parameter curve while Figure 4 represents the ‘Start to middle’ (10–80 N of FF) LBI parameter curve. As a result, ‘Zesy002’ data may lack the 1^st phase of softening data and hence the models reflect more on the exponential FF decrease during the rapid 2^nd phase of softening. For ‘Hayward’ kiwifruit, the data collection may lack the 3rd phase of softening to obtain a wider range of the ripening process. Alternatively, together Figures 3 and 4, a 3-phase curve with a lag phase, an increasing phase and a lag phase is expected and corresponds to the softening curve.

Although in practice kiwifruit ripening is likely to correspond to only part of the LBI parameter time-series curve, a complete time-series curve is desired for developing the calibration model further. Therefore, to gain a better forecasting result of LBI parameters during ripening, fruit at a firmer range with a longer storage time are required to be studied. For time series analysis, more data points (observations) could enhance the model's performance. More data points could be achieved with a shorter data collection interval using an automated system or through a longer storage time using kiwifruit samples during cool storage. However, the model provides an approach to follow fruit in shelf life to provide some decision support for retailers and consumers.

Preliminary research was conducted for collecting kiwifruit LBI parameter data at higher data point resolution by repeated measurement of the same spot at 1-hour intervals. The signal to noise was found to be minimal from the attenuation profile considering the signal–noise ratio of the camera is 37.84 dB. However, fruit skin started to exhibit black spots due to the thermal effect of the laser energy after 60 measurements each at 5 wavelengths.

For within fruit variation, the difference in LBI parameters measured from the shoulder and the equator of kiwifruit was observed due to the curvature of the fruit. In this study, the measurement deviation due to within-fruit variation was minimised by measuring the same imaging region with marked labels and guidance from the pilot laser centred at the labels.

On the other hand, other fruit attributes apart from FF may influence time series analysis of ‘Hayward’ softening. The hairy surface of ‘Hayward’ and the subsurface structure are also likely to contribute to a higher scattering, resulting in lower LBI parameters. In general, ‘Hayward’ has a less uniform shape than ‘Zesy002’. The time series model may be influenced by kiwifruit LBI data obtained from various curvatures at the equatorial region being averaged since a quantitative analysis of the curvature effect was not taken into consideration in the calculation of LBI parameters. The curvature effect from the shape of the fruit was not taken into consideration (Qing et al. 2007; Qin and Lu 2008). Peng and Lu (2006) calculated a sample size correction factor $\sin \theta$ using spatially resolved diffuse reflectance images and found the corrected µ_s’ would be less comparable with the estimated value because $\sin \theta$ < 1, thus the actual µ_s’ should be lower than the estimated µ_s’ at the surface with a larger curvature. In this study, ‘Zesy002’ kiwifruit were graded by industrial export standards and, therefore, are considered to have greater consistency in shape and curvature. However, this was not the case for ‘Hayward’, thus the forecasting error in LBI of ‘Hayward’ may be attributed to both the ripening process and a shape effect.

Moreover, each kiwifruit sample may differ in DMC, SSC and skin colour as they had different pre-harvest conditions (Hewett et al. 1999). Pre-harvest conditions include (micro) climate, irrigation, nutrition, fruit position and orchard management amongst others (Hewett et al. 1999). Differences in LBI parameters among GLs could also be expected when sourced from different GLs despite the kiwifruit being selected with similar FF and stored at the same condition. For ‘Zesy002’ kiwifruit, the trendlines of GLs 5 and 6 appeared largely separated. Similarly, for ‘Hayward’, GLs 1 and 10 had lower LBI parameter values compared with GLs 2 and 8. Therefore, the time series analysis model may be influenced when averaging different trends between growers. However, the sample size of each GL in this study is small due to the measurement capability. In the future, multi-time series analysis may be studied. This allows for analysis of single GL LBI trends and comparisons among GLs.

Conclusion

Kiwifruit LBI parameters through ripening were analysed using time series analysis and increases in LBI parameters were observed. Among the three time series analysis models tested, the α and β coefficients in HWS model are potentially useful for the LBI parameter forecast. Time-series analysis may indicate changes and trends of LBI parameters (DIP, FWHM, SLP and HWQM) which may be associated with kiwifruit ripening. In addition, coefficients in time series analysis models showed the potential to describe the LBI increase rate resulting from different kiwifruit ripeness patterns. However, correlations between LBI and shelf life as well as time series analysis coefficients require further study. As a non-destructive method, LBI parameters and subsequent time-series analyses could be utilised to monitor the kiwifruit ripening progress.

Acknowledgements

The authors would like to acknowledge Christian Regen from ATB, Ian Thomas, John Pedley and Peter Jeffery from Massey University for supporting initial LBI system set-up and system upgrading; Sue Nicholson from Massey University for coordinating after-hour postharvest lab access; Steve Glasgow, Petru Daniels and Evelyn Lupton from Massey University for assisting microscopy sample preparation; Peter Jeffery from Massey University, Zespri Innovation Team and Anthony Pangborn from ARM for organising kiwifruit.

Disclosure statement

No potential conflict of interest was reported by the authors.