An empirical model to estimate the growth of ice crystals for storage of Tilapia at variable temperature conditions

According to Kaale & Eikevik (2013) recrystallization consists of the fusion of small crystals and subsequent solidification in others with the growth of large crystals. This phenomenon affects the quality of frozen products since small crystals are desirable and large crystals generally produce structural cellular damage during frozen storage. Some types of recrystallization naturally occur at constant temperatures because water vapor tends to transfer from high vapor pressure regions. Nevertheless, most problems related to the recrystallization process are the result of temperature fluctuations (Petzold & Aguilera, 2009). Temperature fluctuations during storage have been pointed out by several authors as the main cause of recrystallization in both animal and plant products (Reno et al., 2011).


Introduction
According to Kaale & Eikevik (2013) recrystallization consists of the fusion of small crystals and subsequent solidification in others with the growth of large crystals. This phenomenon affects the quality of frozen products since small crystals are desirable and large crystals generally produce structural cellular damage during frozen storage. Some types of recrystallization naturally occur at constant temperatures because water vapor tends to transfer from high vapor pressure regions. Nevertheless, most problems related to the recrystallization process are the result of temperature fluctuations (Petzold & Aguilera, 2009). Temperature fluctuations during storage have been pointed out by several authors as the main cause of recrystallization in both animal and plant products (Reno et al., 2011). Hartel (1998) defines recrystallization as the melting of smaller sized ice crystals as a separate category, since it can also cause significant changes in the size distribution of ice crystals. This type of recrystallization occurs specifically due to temperature oscillations and involves increasing and decreasing size due to partial melting of the crystal. This condition must occur to a greater extent at higher temperatures and more rapidly for smaller crystals.
On this sense, other authors such as Cooper et al. (2008) have developed a simple unifying model for the crystallization and melting temperatures, showing that homogeneous nucleation and phase transformations driven by the thickening of preexisting surface layers are limiting conditions during this process. Since temperature variation during storage was identified by other authors as a condition that enhances the recrystallization process, it is important to mention that as in most chemical reactions, temperature, activation energy and process constant (recrystallization) are variables that can be related using known mathematical models, and also using more complex numerical methods.
The use of empirical models to represent the recrystallization process in stored frozen food is a simple strategy and limited by the process complexity, but it has the advantage of summarizing a large amount of data in terms of few parameters besides being easier to be constructed. However, since empirical models take the form of the relation function between input and output variables, usually, there is no theoretical basis for this relationship. The empirical model has better adjustment when used within the range of experimental values from which they were generated (Özilgen, 2011).
According to Almonacid-Merino & Torres (1993), fluctuating temperature conditions require continuous temperature measurements in food or heat transfer models to estimate food temperature as a function of time and location. Additionally, ambient temperature records and the thermal properties of foods, as well as the packaging, are required for these estimates. However, Kumar & Panigrahi (2009) argue that the main difficulty in the numerical solution of the heat transport equation is to deal with the great latent heat, which evolves over a very small temperature range. Zuritz & Singh (1985) argue that numerical techniques are necessary to solve thermal energy equations involving phase change for finite objects with temperature-dependent physical properties and non-linear contour conditions. Recent studies such as that conducted by Margeirsson et al. (2012) have been able to develop and validate a 3-D heat transfer model that can be used to predict temperature changes of frozen cod packed in EPS boxes. The numerical model developed can be used to improve the insulation package design or to estimate the shelf life in a profitable way, since these concepts are highly related to the product temperature and can be related to a model of shelf life prediction.
All authors agree that in order to minimize quality losses during processing/storage and to predict shelf life, quantitative kinetic models are needed to express the functional relationship between composition and environmental factors in food quality. The applicability of these models is based on the accuracy of the model and its parameters.
According to Baik et al. (1997), the recrystallization process in vegetable starch can be analyzed by the Johnsson-Mehl-Avrami-Kolmogorov (JMAK) kinetic model according to Equation 1: where: θ is the recrystallized part at time t, E o is the enthalpy at time 0, E t is the enthalpy at time t, E L is the maximum enthalpy, k is the constant referring to the rate at which the process [time -1 ] occurs and n is the Avrami exponent.
The influence of temperature on the starch recrystallization rate can be represented by Q 10 as shown in Equation 2: where: T θ is the rate of recrystallization process depending on the temperature in K.
Additionally, the dependence of the recrystallization process with temperature can be represented by the Arrhenius model shown in Equation 3: where: k is the kinetic constant [T -1 ], A is the frequency factor, E a is the activation energy [cal/mol], t is the temperature [K] and R is the universal gas constant [1.987 cal/mol K].
Based on the activation energy value (E a ), one can derive information about a reaction mechanism or, in some cases, extrapolate data to new conditions; however, extrapolation or interpolation to other temperature values must be performed carefully, respecting conditions and environmental factors in which the process occurs.
Based on the above-mentioned models, the Arrhenius model has the advantage of allowing the extrapolation of data corresponding to a given reaction rate as a function of temperature as long as the activation energy of this process is known. Thus, this model can be used when the degradation process is predetermined by a chemical reaction that determines the reaction rate Petrou et al. (2002). With this background and considering the kinetics of ice crystals growth, it is possible to affirm that the Arrhenius model can be applied to estimate the result of a recrystallization process as a function of different storage temperature conditions. On general, the loss of quality of a food or shelf life can be evaluated by measuring a quality index called A. The variation can usually be represented by the kinetic model: where: k is the rate constant dependent on temperature, product and package characteristics; n is a power factor called reaction order and defines whether the variation rate depends on the amount of A. Of the environmental factors are kept constant, n also determines the shape of the deterioration curve Labuza (1984).
According to Erickson & Hung (1997), the Equation 4 can be written as: where: F (A) is the quality function and k and t are the same parameters of Equation 4. The shape of F (A) depends on the n value. When n is zero, it is called zero-order reaction kinetics, which implies that the quality loss rate is constant under constant environmental conditions. Of n is equal to 1, it is called first-order reaction kinetics, which results in an exponential decrease in the loss rate as the quality decreases. Equation 5 describes the deterioration of any quality attribute as a function of a temperature. This mathematical model obeys to a kinetic of chemical reaction, which can be modeled from empirical data.

Statistical analyses
For this work, the measurements of frozen Tilapia muscle micrographs stored at -18 °C under 3 levels of temperature oscillations (± 3, ± 5, ± 7 and ± 0 °C) published by Gutiérrez et al. (2017) were used. Equivalent diameter histories measured during 90 days of storage were used as a database for the elaboration of the proposed model. The observation of the microstructure size (pores) resulting from temperature fluctuation treatments in relation to the number of days of storage was verified using one-way ANDVA with the aid of Minitab  software 17. The equivalent diameter values for each sampling week and for each of the treatments were expressed with the resulting mean value and corresponding standard deviation. With this statistical treatment, we tried to verify if there were significant differences among the oscillation levels proposed, and then, the Tukey test was used to identify which treatments differed from each other at statistical significance level of p <0.05.

Determination of the recrystallization constant
From the model proposed by Labuza (1984) to predict changes in the properties of a food (Equation 4) as a function of time and environmental conditions, considering the power factor of the reaction order n = 1, the variation of property A is represented by an exponential variation that directly depends on the amount of A. Thus, the following expression can be deduced to evaluate the growth of ice crystals as a function of the recrystallization constant: Ontegrating, we have where: o D : Onitial diameter of ice crystals

Statistical treatment
Since ANDVA showed that there are differences between the corresponding means, the Tukey test was applied to determine the difference between each of the equivalent diameter results obtained from each oscillation level. The test results can be seen in Table 1.

Recrystallization constant
Equation 8 represents the estimated model for the recrystallization kinetics of experimental data. The adjustment of experimental data to the chemical reaction kinetics led to the identification that the growth of ice crystals during the recrystallization process follows a first-order model. This model allowed calculating the kinetic constant (recrystallization constant) k for each of the proposed treatments. Experimental values of the initial mean equivalent diameter (D 0 ) and mean equivalent diameter (Deq) values measured up to the 18 th day of storage were considered for the determination of the recrystallization constant. This period corresponds to the effective time where an increase in the crystal diameter was observed (recrystallization process).
A summary of recrystallization constants obtained is shown in Table 2.
Compared with the results obtained by Monzon Davila (2014) for ultra-fast freezing of gelatin gel model solutions, the k values obtained in this work were very close. These authors report k values of 0.0132 and 0.0242 μm2/day for samples stored at -18 ± 2.5 ° C and -12 ± 2.5 °C, respectively.

Determination of the activation energy (Ea) for the recrystallization process
The activation energy E a corresponding to the recrystallization process can be estimated using the k values obtained in Table 2 With parameters obtained from the adjustment presented in Figure 1, the activation energy for the recrystallization process of Tilapia muscle (Oreochromis niloticus) stored under different temperature oscillation levels was calculated, and the value obtained was: 4758.737 cal/mol.

Application of the Arrhenius model for predicting the growth of ice crystals (Deq) as a function of temperature oscillations
Since the recrystallization process is highly influenced by temperature oscillations during the storage period; and that the Arrhenius equation is based on a relationship between The linearized Arrhenius equation for two pairs of temperature values (T 1 and T 2 ) and the kinetic constant (k 1 and k 2 ) can be used to solve problems of constant prediction as a function of temperature, when one of the points of the straight line is known.
Thus, considering that in the Arrhenius equation only the kinetic constant (k) is temperature dependent, the k 1 and k 2 of Equation 9 can be replaced by Deq 1 and Deq 2 as shown in Equation 10.  Equation 10 represents a prediction of the growth of ice crystals in Tilapia muscle stored under the conditions presented in this work, so it should be emphasized that the validity of this empirical model is conditioned by the following contour conditions: • Non-isothermal conditions during the recrystallization process; • Temperature T 1 = -18 °C (255.15 K) and the mean equivalent diameter measured at that temperature (Deq 1 ) are used as references in the calculation and correspond to the lowest oscillation condition. On other words, the model compares an ideal storage condition and from this, it estimates the expected growth of crystals for a specific day (Deq 2 ), as a function of a less efficient storage condition (T 2 ); • The temperature representative of the oscillation level should be chosen among values of -15, -13 and -11 °C (258.5, 260.5 and 262.15 K), since the experimental data used to generate the model were measured at these temperatures. Dther oscillation temperatures can be used provided they are within the range mentioned above; • The absolute temperature scale must be used for calculations; • The activation energy value should be the experimentally estimated value for this process: 4758.737 cal/mol. Figure 2 shows the equivalent diameter values obtained with Equation 10 for 80 days of storage with temperature oscillations corresponding to ± 3, ± 5 and ± 7 °C (i.e., temperatures representative of oscillation levels -15, -13 and -11 °C) and mean Deq values experimentally obtained under the same oscillation conditions. For the calculation, the mean equivalent diameter value (Deq1) corresponding to treatment Td (storage temperature: -18 ± 0.0 °C) was used as reference.
Since the Arrhenius equation is based on a relationship between reaction rate and temperature, using Equation 10 to estimate the growth of crystals as a function of temperature seems obvious; however, it should be mentioned that the validity of this model is from the contour conditions considered for its formulation (Fernández et al., 2002). These conditions are necessary since the Arrhenius model individually relates the chemical reaction rate with point temperatures; thus, it was considered that, when the temperature oscillation occurred, two conditions were created, one of upper oscillation and one of lower oscillation (temperature representative of the degenerative process). Thus, considering that the recrystallization process is enhanced by differences between freezing temperatures and storage temperatures, it is known that damages at intracellular level occur by the fusion of crystals in the presence of high temperatures; therefore, it could be inferred that the degenerative process will begin at these temperatures.
Predicting the theoretical value of the equivalent diameter of ice crystals during storage is very important since it would allow to estimating changes in product quality and even in shelf life, even in an indirect way. On this perspective, Equation 10 represents a specific model for the specific contour conditions and proposed controlled temperature oscillations. On theory, the proposed model would be valid for any temperature value, since the Arrhenius model is valid in most temperature-dependent chemical kinetics reactions; however, several authors have reported difficulty in modeling the effect of temperature oscillations due to the many secondary reactions that may occur. Nevertheless, the use of the model proposed in Equation 10 may be sufficiently valid to predict the size of crystals, provided that an accurate selection of the experimental conditions is performed. The model proposed in Equation 10 has the limitation of estimating the mean equivalent diameter for a point storage condition and for each instant. That is, the model does not contemplate the cumulative effect of time on crystal growth. On this context, an additional model was proposed, which is capable of describing the interaction of time and the effect of temperature oscillations on the growth of ice crystals. Since the response variable (mean equivalent diameter Deq) is influenced by more than one factor, the response surface methodology was used for the modeling of experimental data obtained by Equation 10. The Deq values plotted in Figure 2 represent a reliable estimate of the growth of ice crystals based on the empirical determination of the recrystallization constant characteristic of the process; therefore the Deq values estimated by Equation 10 for each storage temperature may be used as a database for modeling crystal growth as a function of temperature oscillations and time. Figure 3 shows the response surface for the growth of the mean equivalent diameter of ice crystals during 80 days of storage and storage temperature of -18 °C with 4 temperature oscillation levels ± 0, ± 3, ± 5 and ± 7 °C. Figure 3 shows the influence of temperature oscillations on the growth of ice crystals during the storage time. A directly proportional relation between the degree of variation of the temperature and the result of increasing crystals is evidenced. During the 80 days of storage this behavior is maintained at all levels of oscillation. High levels of growth of the mean equivalent diameter of ice crystals are observed at the storage temperature corresponding to + 7 °C.
With the aid of the Drigin  8.0 software, it was possible to perform the non-linear adjustment of the surface to obtain the corresponding model. The model selected for the surface adjustment was the 2D parabola type as shown in Equation 11.  Table 3.
The satisfactory adjustment of the model and the convergence of data were verified by the correlation coefficient R 2 = 0.917. Thus, it could be concluded that the good adjustment obtained by the model make it a useful tool to describe the recrystallization process in a transient way and resulting from poor storage; provided that the storage temperature has been programmed at -18 °C and oscillating at levels of ± 3, ± 5 and ± 7 °C. However, it should be emphasized that obtaining this model responds to specific experimental conditions; therefore, it is not possible to extend its use to any storage condition.

Conclusions
Predicting the theoretical value of the equivalent diameter of ice crystals during storage is very important, as it enables estimating changes in product quality and even in shelf life, even in an indirect way. On this perspective, the proposed model represents a specific model for specific contour conditions and proposed controlled temperature oscillations. On theory, the proposed model would be valid for any temperature value, since the Arrhenius model is valid in most temperature-dependent   chemical kinetics reactions; however, several authors have reported difficulty in modeling the effect of temperature oscillations due to the many secondary reactions that may occur. However, the empirical model developed has the advantage of being simple and was adequate to predict the growth of ice crystals transiently in frozen Tilapia fillet samples stored at -18 °C, oscillating at levels of ± 0, ± 3, ± 5 and ± 7 °C. The good adjust of model R 2 = 0.917 makes it a useful tool to describe the recrystallization process, as long as the selection of storage temperature and oscillation conditions is carefully performed.