foods Review Rennet-Induced Casein Micelle Aggregation Models: A Review Daniel Salvador 1 , Yoseli Acosta 2 , Anna Zamora 3 and Manuel Castillo 3, * 1 Department of Agroindustrial Science, National University of Trujillo, AV. Juan Pablo II s/n—University City, Trujillo 13011, Peru; dsalvador@unitru.edu.pe 2 School of Agroindustrial Engineering, National University of Trujillo, AV. Juan Pablo II s/n—University City, Trujillo 13011, Peru; yacosta@unitru.edu.pe 3 Department of Animal and Food Science, Centre d’Innovació, Recerca i Transferència en Tecnologia dels Aliments (CIRTTA), Universitat Autònoma de Barcelona, Travessera dels Turons s/n, Bellaterra, 08193 Barcelona, Spain; anna.zamora@uab.cat * Correspondence: manuel.castillo@uab.es; Tel.: +34-93-581-1123 Abstract: Two phases are generally recognized in the enzymatic coagulation of milk: hydrolysis and aggregation, although nowadays more and more researchers consider the non-enzymatic phase to actually be a stage of gel formation made up of two sub-stages: micellar aggregation and hardening of the three-dimensional network of para-κ-casein. To evaluate this controversy, the main descriptive models have been reviewed. Most of them can only model micellar aggregation, without modeling the hardening stage. Some are not generalizable enough. However, more recent models have been proposed, applicable to a wide range of conditions, which could differentiate both substages. Manu- facturing quality enzymatic cheeses in a cost-effective and consistent manner requires effective control of coagulation, which implies studying the non-enzymatic sub-stages of coagulation separately, as numerous studies require specific measurement methods for each of them. Some authors have  recently reviewed the micellar aggregation models, but without differentiating it from hardening.  Therefore, a review of the proposed models is necessary, as coagulation cannot be controlled without Citation: Salvador, D.; Acosta, Y.; knowing its mechanisms and the stages that constitute it. Zamora, A.; Castillo, M. Rennet-Induced Casein Micelle Keywords: enzymatic coagulation; hydrolysis; κ-casein; clotting time; gelation time; coagulation Aggregation Models: A Review. time; micellar aggregation; gel assembly; gel firming; gel hardening Foods 2022, 11, 1243. https:// doi.org/10.3390/foods11091243 Academic Editors: Maria de Lurdes Dapkevicius and Alfredo E. S. 1. Introduction de Borba The coagulation or gelation of milk after the addition of a proteolytic enzyme (e.g., Received: 20 March 2022 chymosin) is usually divided into two overlapping phases: the enzymatic phase and the Accepted: 12 April 2022 non-enzymatic or aggregation phase [1]. This division into two stages has been used by Published: 26 April 2022 various authors [2–5]. In the first stage, the hydrolysis of κ-casein, located on the micellar Publisher’s Note: MDPI stays neutral surface, destabilizes the micelles both electrostatically and sterically, triggering, in the with regard to jurisdictional claims in presence of Ca2+ , the micellar aggregation process [6]. Hardening is typically considered an published maps and institutional affil- extension of aggregation, but, indeed, it seems not to be the case, as the study of this step iations. actually requires a different approach and instrumentation compared to aggregation. In fact, gel hardening is the most complex and currently least understood aspect of enzymatic milk coagulation. Various instruments are available to study gel rheological properties that are critical for cheese manufacturing control, such as the rheologically determined gelation Copyright: © 2022 by the authors. time or the appropriate elastic modulus for cutting [7]. As a result, there is particular Licensee MDPI, Basel, Switzerland. interest in developing methods that can be used to study/quantify gelation in the cheese This article is an open access article vat [8], as control of the gel hardness is critical for the objective selection of the optimal curd distributed under the terms and cutting time that affects both cheese quality and yield. Figure 1 outlines the coagulation conditions of the Creative Commons process of milk. Attribution (CC BY) license (https:// Payne and Castillo [9], who studied milk coagulation using an NIR light scattering creativecommons.org/licenses/by/ sensor, described in the sigmoidal curve generated by the instrument response two periods 4.0/). Foods 2022, 11, 1243. https://doi.org/10.3390/foods11091243 https://www.mdpi.com/journal/foods Foods 2022, 11, 1243 2 of 15 Foods 2022, 11, x FOR PEER REVIEW 2 of 15 that overlap during coagulation (Figure 2): (A) the period corresponding to the phase of enzymatic hydrolysis of micellar destabilization, which goes from the addition of the enzyme to approximately the inflection point of the light scattering ratio (R) or time corresponding to the maximum of the first derivative of R (tmax ); and (B) the period corresponding to the non-enzymatic phase of coagulation which consists of two phases: (a) the period corresponding to the aggregation phase of the destabilized micelles that goes from tmax until the moment in which the gel begins to offer mechanical resistance (e.g., the point that could be identified, approximately, as the time at which a value of G0 = 1 Pa is reached, determined by a rheometer [10,11], related to t2min2 [12]) and (b) the asymptotic period corresponding to the hardening phase of the casein gel that continues until the optimum curd cutting point is reached. Harboe et al. [13] and Fox et al. [8] also describe differences between the aggregation phase and the hardening phase, but from the point of view of changes in viscosity and rheological parameters. Horne and Lucey [14] carried out the last review of the aggregation models but, as mentioned above, without differentiating the aggregation phase from the hardening phase. As seen above, the controversy raised by the separation of the non-enzymatic phase of milk coagulation into phases of aggregation Foods 2022, 11, x FOR PEER REVIEW 2 of 15 and hardening Figure will serve 1. Milk coagulation process.as the basisfrom Reproduced for [5]. studying, in a first approach, the different models that describe the phenomenon of micellar aggregation. Payne and Castillo [9], who studied milk coagulation using an NIR light scattering sensor, described in the sigmoidal curve generated by the instrument response two peri- ods that overlap during coagulation (Figure 2): (A) the period corresponding to the phase of enzymatic hydrolysis of micellar destabilization, which goes from the addition of the enzyme to approximately the inflection point of the light scattering ratio (R) or time cor- responding to the maximum of the first derivative of R (tmax); and (B) the period corre- sponding to the non-enzymatic phase of coagulation which consists of two phases: (a) the period corresponding to the aggregation phase of the destabilized micelles that goes from tmax until the moment in which the gel begins to offer mechanical resistance (e.g., the point that could be identified, approximately, as the time at which a value of G’ = 1 Pa is reached, determined by a rheometer [10,11], related to t2min2 [12]) and (b) the asymptotic period corresponding to the hardening phase of the casein gel that continues until the optimum curd cutting point is reached. Harboe et al. [13] and Fox et al. [8] also describe differences between the aggregation phase and the hardening phase, but from the point of view of changes in viscosity and rheological parameters. Horne and Lucey [14] carried out the last review of the aggregation models but, as mentioned above, without differentiating the aggregation phase from the hardening phase. As seen above, the controversy raised by the separation of the non-enzymatic phase of milk coagulation into phases of aggregation Figure 1. Milk Figure and hardeningcoagulation 1. Milk will process. coagulation serve Reproduced as theprocess. from [5]. in awith basis forReproduced studying, permission first approach, thefrom [Lucey, different J.], [Rennet-Induced models that describe the Coagulation of phenomenon of micellarofaggregation. Milk. In Encyclo-pedia Dairy Sciences]; published by [Elsevier], [2011]. Payne and Castillo [9], who studied milk coagulation using an NIR light scattering sensor, described in the sigmoidal curve generated by the instrument response two peri- ods that overlap during coagulation (Figure 2): (A) the period corresponding to the phase of enzymatic hydrolysis of micellar destabilization, which goes from the addition of the enzyme to approximately the inflection point of the light scattering ratio (R) or time cor- responding to the maximum of the first derivative of R (tmax); and (B) the period corre- sponding to the non-enzymatic phase of coagulation which consists of two phases: (a) the period corresponding to the aggregation phase of the destabilized micelles that goes from tmax until the moment in which the gel begins to offer mechanical resistance (e.g., the point that could be identified, approximately, as the time at which a value of G’ = 1 Pa is reached, determined by a rheometer [10,11], related to t2min2 [12]) and (b) the asymptotic period corresponding to the hardening phase of the casein gel that continues until the optimum curd cutting point is reached. Harboe et al. [13] and Fox et al. [8] also describe differences between the aggregation phase and the hardening phase, but from the point of view of changes in viscosity and rheological parameters. Horne and Lucey [14] carried out the last review of the aggregation models but, as mentioned above, without differentiating the aggregation Figure 2. phase from thecurves Coagulation hardening phase. As seen of pasteurized above, whole the R milk. controversy is the lightraised by scattering ratio, R’ and R” theare separation the firstofand the second non-enzymatic phaseof derivatives ofR,milk coagulationtinto phases respectively. is the of aggregation time to the maximum of the first max and hardening will serve as the basis for studying, in a first approach, the different models derivative of R. t2min2 is the time to the second minimum of the second derivative of R. A is the that describe the phenomenon of micellar aggregation. enzymatic period. B is the nonenzymatic period. a is the aggregation period. b is the hardening period. Reproduced with permission from [Payne, F.A.; Castillo, M.], [Light Backscatter Sensor Applications in Milk Coagulation. In Encyclopedia of Agricultural, Food, and Biological Engineering]; published by [Taylor & Francis], [2007]. Foods 2022, 11, 1243 3 of 15 2. Micellar Aggregation After the enzymatic hydrolysis of κ-casein, the curd formation phase begins, consisting, in the first instance, in the aggregation of hydrolyzed casein micelles. Aggregation of destabilized micelle begins when there is sufficient hydrolyzed κ-casein 65–90%, according to [15], and if the activity of calcium ions and temperature are adequate [16]. Initially, small aggregates are formed in the form of chains (Figure 3), instead of lumps [6]. As more micelles are aggregating, the apparent viscosity and the amount of casein macropeptide (CMP) released are observed to increase. Chain formation could be explained using the geometric model described by Dalgleish and Holt [17], who state that two hydrolyzed casein micelles will only aggregate if they have sufficiently large holes on their surface formed by hydrolyzed κ-casein and if they collide at the proper orientation. The inter-particle bonds could be established through calcium bridges and/or hydrophobic interactions [8]. Modification of histidyl, lysyl, and arginyl residues in κ-casein inhibits the aggregation phase. This finding suggests that a positively charged group on para-κ-casein reacts electrostatically with unidentified negative charged sites. In native micelles, this positive site may be shielded by the CMP segment of κ-casein, becoming exposed and, subsequently, reactive when this peptide is released [18]. Changes in the surface hydrophobicity of casein Foods 2022, 11, x FOR PEER REVIEW micelles during the gelation process have also been confirmed using a fluorescent 4 of 16 marker, 8-anilino naphthalene-1-sulfonate [19,20]. Figure 3. Release of macropeptides, expressed as fraction (•). changes in relative viscosity () during Figure the course 3. Release of rennet coagulation. of macropeptides, expressedOnas top of the fraction figure, casein (●).Changes micelle viscosity in relative aggregation(■)progress dur- ing is the course of rennet graphically coagulation. represented. On top Reproduced of the with figure, casein permission micelle from [Fox, aggregation P.F.; progress Guinee, T.P.; Cogan,isT.M.; graphically McSweeney,represented. Reproduced of P.L.H], [Fundamentals with permission Cheese from Science, 2nd [Chen, ed.]; C.-C.;by published Chen, L.-Y.; Li, [Springer], W.-T.; [2017]. Chang, K.-L.; Kuo, M.-I.; Chen, C.-J.; Hsieh, J.-F], [Influence of Chymosin on Physicochemical and Hydrolysis The Characteristics aggregationofisCasein due inMicelles part toand theIndividual Caseins.]; van der Waals published attraction, by it[Nanomateri- but also requires als],participation [2021]. 2+ of Ca , as the van der Waals attraction is not sufficient. The need of a sufficient amount of Ca2+ in the serum phase of milk clearly shows its role. Presumably, the Figure effect of 4Ca shows the 2+ ions is percentage of CMP double. First, released aselectrostatic Ca2+ decreases a function of time at different repulsion pH by neutralizing values. negative charges on the micelles [21], being more efficient than H ions in the pHthe It can be observed that, as the pH decreases, less time is needed + to release range amount underofconsideration. CMP necessary to cause Second, theions Ca2+ aggregation ofbetween can bridge the micelles [23], sites negative as acidification on paracasein progressively micelles bysolubilizes means of the salt micellar bridges.calcium phosphate, It is very well-knownreleasing ionic calcium, that lowering the pH bring- of milk ing forward the start of the aggregation phase [24]. The Berridge clotting time or visual coagulation time is the time at which the first flocs of hydrolyzed casein are visible on the walls of a glass tube [25,26]. The primary enzyme reaction is ~85% complete at a time equivalent to 60% of the visual clotting time. Between 60 and 80% of the visual coagulation time, the rennet-destabilized micelles start to aggregate gradually [8]. Foods 2022, 11, x FOR PEER REVIEW 4 of 15 Foods 2022, 11, 1243 4 of 15 Figure 3. (a) Release of macropeptides, expressed as fraction (●) and (b) changes in relative viscosity (■) during the course of rennet coagulation. On top of the figure, casein micelle aggregation pro- gress is graphically considerably increases Ca2+ activity represented. Adapted[22] from by[6]. shifting the existing ionic balance between the soluble and micellar phases, releasing H+ from the micelle. Figure Figure44shows showsthethepercentage percentageofofCMPCMPreleased releasedas asaafunction functionofoftime timeatatdifferent differentpHpH values. values.ItItcan canbe beobserved observedthat, that,asasthe thepHpHdecreases, decreases,less lesstime timeisisneeded neededtotorelease releasethe the amount amountofofCMPCMPnecessary necessarytotocause causethetheaggregation aggregationofofthe themicelles micelles[23], [23],as asacidification acidification progressively progressively solubilizes solubilizes the micellar calcium calcium phosphate, phosphate,releasing releasingionic ioniccalcium, calcium,bringing bring- ing forward the start of the aggregation phase forward the start of the aggregation phase [24]. The [24]. The Berridge clotting time time or visual or visual coagulation coagulation time is the time at which the first flocs of hydrolyzed casein visible time is the time at which the first flocs of hydrolyzed casein are on the are visible on walls of a of the walls glass tubetube a glass [25,26]. TheThe [25,26]. primary primaryenzyme enzyme reaction is is reaction ~85% ~85% complete completeatataatime time equivalent equivalentto to60% 60%ofofthe thevisual visualclotting clottingtime. time.Between Between60 60and and80% 80%ofofthe thevisual visualcoagulation coagulation time, time,the therennet-destabilized rennet-destabilizedmicelles micellesstart startto toaggregate aggregategradually gradually[8]. [8]. Releaseofofcasein Figure4.4.Release Figure caseinmacropeptide macropeptideduring duringcoagulation coagulationofofmilk milkatatdifferent differentpH pHvalues. values.Repro- Repro- ducedfrom duced with [23]. permission from [Corredig, M.; Salvatore, E.], [Enzymatic Coagulation of Milk. Advanced Dairy Chemistry: Volume 1B: Proteins: Applied Aspects]; published by [Springer], [2016]. Castillo et al. [27], coagulating goat milk, found that the start of the aggregation phase couldCastillo et al. [27], be objectively coagulating determined goatan using milk,NIR found lightthat the startsensor, scattering of the aggregation phase by the inflection could be objectively determined using a NIR light scattering sensor, by point of the scattering ratio, R, (Figure 2), which corresponded to 88% of the Berridge the inflection point of the scattering ratio, R, (Figure 2), which corresponded to 88% of the coagulation time. Tabayehnejad et al. [28], working with low heat skimmed milk powder Berridge coagulation time. and Tabayehnejad with et al. [28], working the same methodology, with found that thelow heat start skimmed milk of aggregation powder occurs at 84%and ofwith the the same methodology, found that the start of aggregation occurs at Berridge coagulation time, confirming two undeniable facts: that the detection of an ag-84% of the Berridge coagulation gregation time, point confirming mainly depends twoonundeniable the sensitivityfacts: ofthat the detection of the measurement an aggregation technique to the presence of the aggregates and also that the two processes (proteolysis and subsequentof point mainly depends on the sensitivity of the measurement technique to the presence the aggregates aggregation) and also overlap that[4]. in time theTherefore, two processes (proteolysis the value and subsequent of tmax (Figure aggregation) 2) was proposed as an overlap in time [4]. Therefore, the value of t max approximate target value of aggregation start time [28,29]. (Figure 2) was proposed as an approximate target value of aggregation start time [28,29]. 3. Aggregation Models of Hydrolyzed Casein Micelles 3. Aggregation Models of Hydrolyzed Casein Micelles Numerous Numerousefforts efforts have been made have been made to tomathematically mathematicallydescribe describethe theaggregation aggregation reac- reaction tion of hydrolyzed casein micelles. They have been proposed from simpler models, of hydrolyzed casein micelles. They have been proposed from simpler models, funda- fun- damentally empirical,totomore mentally empirical, moresophisticated sophisticated ones, ones, based based on on computational computational simulations simulationsofof Brownian dynamics, which constitute clear evidence of the interest in knowing Brownian dynamics, which constitute clear evidence of the interest in knowing the thebehavior behav- ior of this phase of coagulation. Next, the main models that have been studied/proposed of this phase of coagulation. Next, the main models that have been studied/proposed will will be described. be described. Foods 2022, 11, x FOR PEER REVIEW 5 of 15 Foods 2022, 11, 1243 5 of 15 3.1. Holter Model Almost 150 years ago, Storch and Segelcke (1874) (cited by McMahon et al. [30]) 3.1. Holter Model found that the coagulation time is inversely proportional to the enzyme concentration, accordingAlmost 150formula: to the years ago, Storch and Segelcke (1874) (cited by McMahon et al. [30]) found that the coagulation time is inversely proportional to the enzyme concentration, according to the formula: t E=k (1) tc E = k (1) where tc is the coagulation time, E is the enzyme concentration, and k is the constant. where It was tc is soon found that time, the coagulation this relationship E is the enzyme only concentration, holds in a small andrange k is of thepH, tempera- constant. ture, andIt wasenzyme concentration. soon found For this reason, that this relationship onlyHolter holds in [31] separately a small rangeconsidered the du- of pH, temperature, ration of the hydrolysis and enzyme concentration. phaseFor(ththis ) and the aggregation reason, Holter [31] time (ta), defining separately consideredthe coagulation the duration of the time as hydrolysis phase (th ) and tc = th + ta. Subsequently, the aggregation Foltmann time (taHolter’s [32], assuming ), defining the coagulation proposal, proposed timethe as tc = th equation following + ta . Subsequently, for tc: Foltmann [32], assuming Holter’s proposal, proposed the following equation for tc : C tct == E + C +t (2) ta (2) [E] where wheretataisisthe theaggregation aggregationtime, time,tctis the coagulation time, [E] is the enzyme concentration, c is the coagulation time, [E] is the enzyme concentration, and C is the constant. and C is the constant. InInthis thismodel, model,the theintercept, intercept, ta, twould come to represent the minimum duration of co- a , would come to represent the minimum duration of agulation when the enzyme concentration coagulation when the enzyme concentration tends tends to infinity (i.e., hydrolysis to infinity would be (i.e., hydrolysis instan- would be taneous, instantaneous, not limiting the aggregation rate). As can be seen, this relationship isempir- not limiting the aggregation rate). As can be seen, this relationship is purely purely ical, but it isbut empirical, a relevant relationship it is a relevant that must relationship thatbemust satisfied as a starting be satisfied point for as a starting anyfor point more any descriptive mechanistic models. Figure 5 shows the representation more descriptive mechanistic models. Figure 5 shows the representation of Equation (2). of Equation (2). The equation The equation assumes that there assumes that isthere no overlap between is no overlap proteolysis between and theand proteolysis aggregation phase, the aggregation and that the degree of proteolysis is always the same at the time of phase, and that the degree of proteolysis is always the same at the time of coagulation [33]. coagulation [33]. tc ta 1/E Figure 5. Relationship between enzyme concentration (E), aggregation time (ta ), and coagulation Figure 5. Relationship between enzyme concentration (E), aggregation time (ta), and coagulation time (t ). Reproduced with permission from [Dalgleish, D.G.], [A New Calculation of the Kinetics of time (tc).c Adapted from [34]. the Renneting Reaction.]; published by [J. Dairy Res.], [1988] [34]. 3.2. 3.2.Payens PayensModel Model Payens Payens[35], [35],Payens Payensetetal. al.[36], [36],and andPayens PayensandandWiersma Wiersma[37][37]proposed proposedaamodelmodelforfor aggregation aggregation based on Michaelis–Menten kinetics for hydrolysis and Smoluchowskiaggre- based on Michaelis–Menten kinetics for hydrolysis and Smoluchowski aggre- gation gationkinetics. kinetics. They They considered that the considered that theextent extentofofmicelle micelleaggregation aggregation is is proportional proportional to to the the amount of cleaved κ-casein. Thus, in the early stages of the reaction, small amount of cleaved κ-casein. Thus, in the early stages of the reaction, small concentrations concentra- tions of micelles of micelles that that can aggregate can aggregate maymayhavehave beenbeen originated originated by hydrolysis by hydrolysis [38]. [38]. There- Therefore, fore, therethere willa lag will be be aperiod lag period during during whichwhich no effect no great great of effect the of the enzyme enzyme on the on the average average particle particle size will be observed, proposing the following model for the size will be observed, proposing the following model for the coagulation time: coagulation time: s 2 tc ≈ (3) ks Vmax Foods 2022, 11, 1243 6 of 15 where tc is the coagulation time, ks is the aggregation rate constant, and Vmax is the maximum rate of hydrolysis of κ-casein. This model described well the coagulation kinetics of isolated κ-casein generated with pepsin [39]. The main problem with the model was that the rate constant for aggregation (ks ) of the hydrolyzed micelles was very low compared to the constants found experimentally. This has been attributed to the fact that ks is an average of all values, from very low values at the beginning and high values at the end of the aggregation reaction [3,40–42]. 3.3. Step Function Model Dalgleish [43] proposed the step function model using a binomial distribution. This model is based on the idea that the enzyme would act on the surface of the micelle in a random way and that a critical level of proteolysis must be reached before aggregation is possible. The model assumes that the hydrolysis and aggregation phases can be modeled taking into account that ~97% of the κ-casein in a micelle needs to be hydrolyzed for it to participate in aggregation. Applying a binomial distribution to the hydrolysis phase controlled by Michaelis– Menten kinetics and Smoluchowski aggregation kinetics, we obtain:     kM 1 αc . S0 1 Mcrit tc = th + ta ∼ ln + + −1 (4) Vmax (1 − αc ) Vmax 2ks . C0 M0 where tc is the coagulation time, th is the hydrolysis time, ta is the aggregation time, kM is the Michaelis–Menten constant (it approximates the dissociation constant for the enzyme-substrate complex), Vmax is the maximum rate of proteolysis at infinite substrate concentration, αc is the extent of κ-casein hydrolysis, S0 is the initial concentration of κ- casein, ks : is the aggregation rate constant, C0 is the concentration of aggregating material, Mcrit is the average molecular weight at coagulation time (≈10 micellar units), and M0 is the average molecular weight at t = 0. This calculation method for chymosin-induced casein micelles aggregation is similar to that of Payens [36], but it incorporates certain important differences such as the use of full Michaelis–Menten kinetics and the description of the aggregation mechanisms based on probability functions. The calculations agree well with the experimental evidence, in terms of enzyme variation and substrate concentration [44]. If there is any inconvenience in Dalgleish’s model, it is its complexity, as an analytical solution is not reached, and simplifications are needed to obtain such solutions [3]. 3.4. The Energy Barrier Model The energy barrier model proposed by Darling and van Hooydonk [45] is based on the idea that an energy barrier (Equation (5)) is gradually lowered due to proteolysis, leading to an increasing probability of collision during the reaction of aggregation. The stability factor, W, characterizes the relationship between the particle collision frequency and the frequency with which these collisions are successful, resulting in aggregation. The coagulation of casein micelles by rennet action can then be satisfactorily described by determining the relationship between the stability factor W, which represents a measure of the energy barrier, and the concentration of κ-casein as a function of time [45]. Wt = W0 exp(−Cm .V.t) (5) where Wt is the stability factor at time t for casein micelles, W0 is the initial stability factor for casein micelles, Cm is the constant that relates the stability of the casein micelle to the concentration of κ-casein, V is the enzymatic reaction rate, and t is the time. Foods 2022, 11, 1243 7 of 15 Considering first-order enzyme kinetics (i.e., simplification from Michaelis–Menten kinetics) for the hydrolysis phase and Smoluchowski aggregation kinetics, Equation (6) is obtained for the coagulation time, tc :     1 1 W exp(−Cm .S0 ) 1 1 tc = S0 + (exp(−Cm .S0 ) − 1) + 0 − (6) V Cm ks nc n0 where tc is the coagulation time, S0 is the initial concentration of κ-casein, ks is the aggrega- tion rate constant, nc is the concentration of casein aggregates at time tc , and n0 is the initial concentration of casein micelles. The model describes the variations observed in the coagulation time (tc ) as a function of temperature, rennet concentration, and total protein concentration. Additionally, it demonstrated the existence of a delay time in the start of aggregation equivalent to 60% clotting time [45]. This energy barrier model has been tested extensively [4]. The difficulty with the model lies in its simplification of Michaelis kinetics to allow calculations [39]. 3.5. The Functionality Model Hyslop [46] and Hyslop and Qvist [47] proposed the functionality model, which is based on the fact that the continuous removal of CMP would lead to multiple holes on the surface of the micelle and that, therefore, the aggregation rate constant, kij, should be a function of the number of free reactive sites on the aggregation particles, which met the required conditions for a polyfunctional model of the Flory–Stockmayer type [48]: h i kij = K 4 + 2(f − 2)(i + j) + (f − 2)2 ij (7) where K is the proportionality factor, f is the number of functional sites (functionality), and I and j are the number of aggregating particles of type i and j, respectively. If f = 1, only dimer formation is possible, while if f = 2, the model predicts linear polymers and if f > 2, chain branching occurs and gelation is possible. Initially, at time t = 0, f = 0, the micelles have no reactivity or inclination to aggregate [4]. Although the Stockmayer expression (Equation (7)) does not include explicit geometric assumptions, it is generally assumed that the aggregating entities are large extended polymers in which many functional sites are possible, and in which steric problems are improbable. With spherical particles, such as casein micelles, steric problems could easily appear when two micelles become close enough. Therefore, K, f, or both in Equation (7) would be expected to be low. However, Stockmayer’s expression would still apply [3]. 3.6. Adhesive Hard-Sphere Model Reactions leading to aggregation of hydrolyzed micelles have also been studied using the adhesive hard-sphere theory, proposed by de Kruif [49,50]. Instead of attributing the increase in viscosity to the formation of aggregates, de Kruif proposed a different approach based on considering casein micelles as rigid, sterically stabilized spheres that become sticky, or adhesive, as the κ-casein is being proteolyzed and released. The relative viscosity, nr , of skim milk is described as:   1.9 nr = 1 + 2.5Φ + 5.9 + Φ2 (8) τ where Φ is the particle volume fraction, which varies as hydrolysis and aggregation occurs, and τ is a stickiness parameter related to the depth, ε, of an attractive potential created when κ-casein pili (i.e., CMP) are separated from the casein micelle surface by rennet. Figure 6 shows two interacting casein micelles, where V(r) is the potential of the “square well”, with width ∆ and depth ε, and σ is the diameter of the hard core. Foods Foods 11, x11, 2022, 2022, FOR1243 PEER REVIEW 8 of 815of 15 Figure 6. Two interacting casein micelles. The hard-core repulsion is preceded by an attractive in- teraction modeled as a “square well” potential. The depth of the well depends on the stability of the κ-casein brush. Reproduced from [50]. Chymosin induces an initial decrease in viscosity that exceeds a minimum before an exponential increase takes places as coagulation progresses (Figure 7). Therefore, rota- tional viscometers provide limited information about coagulation, as during this process, the viscosity of milk tends towards infinity at visual coagulation time [51]. The initial decrease in viscosity arises due to the decrease in the particle volume frac- tion, Φ, when the hairy layer of κ-casein is hydrolyzed. This decrease in hydrodynamic size has been measured experimentally using dynamic light scattering techniques in diluted [52] and concentrated [53] micellar suspensions. Subsequently, the viscosity starts to increase due to aggregation and network formation. An important disadvantage of the hard-sphere modelFigure 6. Two interacting casein of micelles. The hard-core repulsion is preceded by an attractive Figureis6. that Twothe time interactingdependence casein micelles. relative viscosity The hard-core is mostly repulsion related is preceded tobythe an proteolysis attractive in- interaction reaction, teraction which modeled modeled limits as a “square its use to as a “square well” preliminary well” potential. potential.stages The depth of of aggregation The depth the well depends [54]. Indeed, of the well depends on the stability on the according stability to of of the the publishedκ-casein κ-casein brush. brush. Reproduced data,Reproduced with from [50]. the hard-sphere permission from [De Kruif, C.G.], [Casein model only seems to describe the coagulation process upMicelle Interactions]; published to the by [Int. Dairytime, visual coagulation J.], [1999]. due to the existence of an asymptote at about 100% visual Chymosin clotting time, ~50induces an initial min in Figure 7. Asdecrease it can bein viscosity observed in that exceedstheory this figure, a minimum before an has successfully exponentialChymosin increase induces takes an initial places as decrease in viscosity coagulation progressesthat(Figure exceeds7). a minimum Therefore, before rota- an estimated viscosity changes up to the point where aggregation occurs [5]. As a result, it can- exponential increase takes places as coagulation progresses (Figure 7). Therefore, rotational tional not viscometers provide informationprovide on the limited information complete about kinetics of coagulation, aggregation and gelas firming, during thistwoprocess, highly viscometers provide limited information about coagulation, as during this process, the the viscosity significant of events milk in tends the towards definition of infinity the gel at visual properties coagulation [4]. time [51]. viscosity of milk tends towards infinity at visual coagulation time [51]. The initial decrease in viscosity arises due to the decrease in the particle volume frac- tion, Φ, when the hairy layer of κ-casein is hydrolyzed. This decrease in hydrodynamic size has been measured experimentally using dynamic light scattering techniques in diluted [52] and concentrated [53] micellar suspensions. Subsequently, the viscosity starts to increase due to aggregation and network formation. An important disadvantage of the hard-sphere model is that the time dependence of relative viscosity is mostly related to the proteolysis reaction, which limits its use to preliminary stages of aggregation [54]. Indeed, according to published data, the hard-sphere model only seems to describe the coagulation process up to the visual coagulation time, due to the existence of an asymptote at about 100% visual clotting time, ~50 min in Figure 7. As it can be observed in this figure, theory has successfully estimated viscosity changes up to the point where aggregation occurs [5]. As a result, it can- not provide information on the complete kinetics of aggregation and gel firming, two highly significant events in the definition of the gel properties [4]. Figure Figure 7. Relative 7. Relative viscosity viscosity of skim of skim milk milk during during the coagulation the coagulation process. process. The drawn The drawn line represents line represents the theoretical the theoretical calculation calculation in which in which thethe micelles micelles areconsidered are consideredtotobe be sticky sticky hard spheres. spheres. Reproduced Repro- duced from [50]. with permission from [De Kruif, C.G.], [Casein Micelle Interactions]; published by [Int. Dairy J.], [1999]. The initial decrease in viscosity arises due to the decrease in the particle volume frac- tion, Φ, when the hairy layer of κ-casein is hydrolyzed. This decrease in hydrodynamic size has been measured experimentally using dynamic light scattering techniques in diluted [52] and concentrated [53] micellar suspensions. Subsequently, the viscosity starts to increase due to aggregation and network formation. An important disadvantage of the hard-sphere model is that the time dependence of relative viscosity is mostly related to the proteolysis reaction, which limits its use to preliminary stages of aggregation [54]. Indeed, according to published data, the hard-sphere model only seems to describe the coagulation process up to the visual coagulation time, due to the existence of an asymptote at about 100% Figure 7. Relative viscosity of skim milk during the coagulation process. The drawn line represents the visual clotting theoretical time, ~50 calculation min in in which theFigure 7. are micelles As considered it can be observed in hard to be sticky this figure, spheres.theory Repro-has successfully duced from [50]. estimated viscosity changes up to the point where aggregation occurs [5]. As Foods 2022, 11, 1243 9 of 15 Foods 2022, 11, x FOR PEER REVIEW a result, it cannot provide information on the complete kinetics of aggregation and gel 9 of 15 firming, two highly significant events in the definition of the gel properties [4]. 3.7. Fractal Aggregation Model 3.7. Fractal BremerAggregation Model et al. [55–57] applied fractal theory to acid-induced casein particle aggregation. It hasBremer also beenet al. [55–57] applied demonstrated thatfractal theory to acid-induced rennet-induced aggregation of casein caseinparticle micellesaggrega- can be tion. It has modeled also beentodemonstrated according that rennet-induced fractal mechanisms aggregation [58–62]. Aggregates can beof casein micelles considered can fractals if be modeled their geometric according to fractal structure mechanisms is similar when viewed [58–62]. overAggregates a reasonably canlarge be considered frac- range of scales (that tals ifis, irrespectively their of the observation geometric structure scale selected). is similar when viewed over As athe number of reasonably particles large range inof an aggregate scales (that is,increase with time, irrespectively of thea observation relationshipscalebetween the size selected). Asof theannumber aggregate and the of particles number in of particles, an aggregate Np , iswith increase obtained time, as follows: a relationship between the size of an aggregate and the number of particles, Np, is obtained as follows:  D R Np = R (9) N =a (9) a where where R is the radius of the aggregate, a is the radius of R is the radius of the aggregate, a is the radius of the the primary primary particle, particle, and and DD is the is the fractal dimension. fractal dimension. Assuming aa sphere Assuming sphere ofof radius radius R R containing containing aa number, number, Ns, Ns, of of closely closely packed packed particles, particles, the volumetric the volumetric fraction fraction of of the the particles particles inin aa fractal fractal aggregate, aggregate, ϕφA , is given by: A, is given by:   D−3 Np N RR ϕA φ= = = = (10) Ns N aa As As DD << 33 invariably invariably (in(in three-dimensional three-dimensional space),space), this this means means that that the the particle particle volume volume fraction fraction ϕφAA will will decrease decrease withwith increasing increasing values values of of R R as as the the group group grows [63]. At grows [63]. At aa certain certain radius (R ), to which the aggregate is known as the fractal ‘blob’, radius (Rg ), to which the aggregate is known as the fractal ‘blob’, ϕA will be equal g φ A will be equal to the to volume fraction of particles in the system φ and all aggregates will the volume fraction of particles in the system ϕ and all aggregates will touch, forming touch, forming a con- a tinuous continuousgel gel network network [5].[5]. Bremer Bremer et et al.al. [55] defined [55] definedthethegelation gelationpoint pointbase baseon onthis this event, event, which infers which infersthatthatall allthe theparticles particlespresent present in in thethe system system areare incorporated incorporated intointo the aggre- the aggregates. gates. The The real real question question is whether is whether this corresponds this corresponds to the experimentally to the experimentally measured measured rhe- rheological ological gel pointgel point (crossover (crossover point—when point—when the elastic the elastic G modulus, 0 , becomes modulus, G′,higher becomes that higher that the viscous the viscousG00modulus, modulus, , which isG″, which isgeneralized sometimes sometimes generalized empirically as G0 = 1 Pa);asthis empirically G′ =seems 1 Pa);tothis be seems to be the assumption made [4]. Furthermore, assuming that the assumption made [4]. Furthermore, assuming that the fractal blob coincides with the the fractal blob coin- cides moment withinthe moment which in which the the rheometer rheometer yields yields arheological a measurable measurableresponse, rheological theresponse, obvious the obvious question is howquestion is howinthe the increase the increase in the elastic elastic modulus modulus takes of the structure of theplaces, structure takes according to the fractal places, accordingmodelto [64]. the fractal model [64]. According to [4], the fractal network structure is similar to the image represented in Figure 8. Formationofofa aweak Figure 8. Formation weak gelgel in in which which most most of the of the particles particles are weakly are weakly bound bound to thetonetwork. the net- LargerLarger work. circles circles represent the fractal represent ‘blob’ concept. the fractal Reproduced ‘blob’ concept. from Horne Reproduced and Bank [4]. with permission from [Horne, D.S.], [Banks, J. Rennet-induced Coagulation of Milk. In Cheese Chemistry, Physics and Microbiology]; Walstra published [63] mentioned by [Elsevier], [2004]. that an intermicellar rearrangement is also possible as the aggregates are formed (short-term rearrangement); the particles wrap around each other until they have obtained a higher coordination number, which indicates a more stable configuration (Figure 9). Foods 2022, 11, 1243 10 of 15 Walstra [63] mentioned that an intermicellar rearrangement is also possible as the Foods 2022, 11, x FOR PEER REVIEW aggregates are formed (short-term rearrangement); the particles wrap around each 10 ofother 15 until they have obtained a higher coordination number, which indicates a more stable configuration (Figure 9). Figure 9. Short-term rearrangement. (a) Particle rearrangement. (b) Fractal group in two dimensions, Figure 9. Short-term rearrangement. (a) Particle rearrangement. (b) Fractal group in two dimen- where a short-term sions, where rearrangement a short-term rearrangementhashas occurred. Reproduced occurred. Reproducedwith permission from [63]. from [Walstra, P.], [Physical Chemistry of Foods]; published by [Marcel Dekker], [2003]. The fractal dimension of casein gels obtained by enzymatic coagulation can be ob- The fractal dimension of casein gels obtained by enzymatic coagulation can be obtained tained using techniques such as permeametry, turbidimetry, or rheology, reporting values using techniques such as permeametry, turbidimetry, or rheology, reporting values for this for this type of gel from 2.17 to 2.4 [44]. type of gel from 2.17 to 2.4 [44]. The fractal model can describe the aggregation of casein micelles until the formation The fractal model can describe the aggregation of casein micelles until the formation of the network occurs, at which point the fractal structure would slowly disappear due to of the network occurs, at which point the fractal structure would slowly disappear due to a a massive restructuring of the aggregates, producing greater compaction of the aggregates massive restructuring of the aggregates, producing greater compaction of the aggregates and creating larger pores in the gel [65,66]. It is therefore not possible to efficiently model and creating larger pores in the gel [65,66]. It is therefore not possible to efficiently model the hardening stage of the gel formation. the hardening stage of the gel formation. The various models mentioned throughout this review essentially describe the in- The various models mentioned throughout this review essentially describe the increase crease in average in average molecular molecular massmass of theofmicellar the micellar aggregates aggregates as a function as a function of The of time. time.average The average molecular mass is easily measurable by static light scattering techniques molecular mass is easily measurable by static light scattering techniques [4]. However, due[4]. How- ever, dueproblems to the to the problems encountered encountered with multiple with multiple light scattering, light scattering, where the where the detected detected photon are photon scattered several times as they travel through the suspension media due to thedue are scattered several times as they travel through the suspension media to the of presence presence multiple scatterers, these techniques are applicable only to highly dilute suspensions sus- of multiple scatterers, these techniques are applicable only to highly dilute [67] or pensions [67] or over very short path lengths, as in the turbidity measurements over very short path lengths, as in the turbidity measurements [68]. A main disadvantage [68]. A main disadvantage of all these models is that most describe the total or partial hydrolysis of all these models is that most describe the total or partial hydrolysis and aggregation of and aggregation casein micelles,ofbut casein in nomicelles, case the but in no case complete the largely process, complete process, due largely related to limitations due to to limitations related to experimental experimental techniques. techniques. 3.8.3.8. Light Scattering Light Model Scattering Model AA method method that hashas that apparently apparently overcome overcome thethe aforementioned aforementioned drawbacks drawbacks is the one is the one proposed proposed byby Castilloetetal.al.[29], Castillo [29],describing describingthe the kinetics kinetics of light scattering scatteringchange changeduring duringthe theaggregation aggregationsteps of of steps hydrolyzed hydrolyzed casein micelles casein andand micelles subsequent gel hardening. subsequent gel hardening.ThisThis model couldcould model differentiate the aggregation differentiate the aggregation stagestage from from the hardening stage.stage. the hardening It combines second It combines orderorder second kinetics to describe kinetics to describe the aggregation the aggregation of casein micelles of casein andand micelles first order first kinetics order that kinetics would that would correspond correspond to to casein casein gelgel hardening. hardening. In In particular, particular, thethe model model considersthe considers thefollowing followingsecond-order second-orderreaction reactionforforthe theaggre- aggrega- tion gation stage: R∞A − R0A   Ra = R∞A − R −R (11) R = R − 1 + (R∞A − R0A )k2 (t0 − tmax ) (11) 1 + (R − R )k (t − t ) where Ra is the diffuse reflectance ratio attributed to the aggregation period; R∞A is where Ra is the the diffuse diffuse reflectance reflectance ratio at t = ratio attributed ∞, during theto the aggregation aggregation period;period; R0A isRthe ∞A is the diffuse diffuse reflectance ratio at t = ∞, during the aggregation period; R0A is the diffuse reflec- tance ratio at tmax; k2 is the second order reaction rate constant; t′ is the time after enzyme addition; and tmax is the starting time for the aggregation reaction. Applying light scattering to monitor milk coagulation by rennet, Castillo et al. [29] observed a systematic distortion of the scattered light profile during coagulation. At low protein concentration and temperature, a distinct shoulder separated the light scattered coagulation profile, which tended to disappear with increasing temperature and protein concentration (Figure 10). This behavior was attributed to the presence of two distinct but Foods 2022, 11, 1243 partially overlapping phenomena: the aggregation of casein micelles and the hardening 11 of 15 of the three-dimensional gel structure. Apparently, these two stages of gel formation were not clearly distinguished at high protein concentration and temperature, possibly because these coagulation reflectance ratio atconditions tmax ; k2 isfacilitated the second theorder overlap of these reaction ratementioned 0 is the time after constant; tphenomena. enzymeThe addition; possible aggregation and tmax is themechanisms involved starting time for thein casein gelsreaction. aggregation are now being studied usingApplying computational simulations light scattering of Brownian to monitor dynamics [5,69]. milk coagulation WithCastillo by rennet, this type of [29] et al. ap- proach, observed the role of repulsive a systematic barriers distortion of theand different scattered lightbonds profiletypes, as coagulation. during well as the behav- At low ior/movements/collisions protein concentration andoftemperature, the particle on the aggregation a distinct shouldermechanisms separated the and kinetics, light can scattered be studied and coagulation better profile, understood. which tended toIndisappear other words, withthis technique increasing is very useful temperature in our and protein understanding of how10). concentration (Figure particle interactions This behavior wasinfluence attributedgelation properties to the presence [5]. distinct of two New tech- but partially niques foroverlapping monitoringphenomena: milk coagulationthe aggregation of casein are also being tested,micelles and thedensity such as photon hardeningwave of the three-dimensional spectroscopy (PDW) [70],gelfast structure. FourierApparently, these two stages transform continuous of gel formation cyclic voltammetry [71],were and not clearly distinguished fluorescence at high protein [72], demonstrating concentration the continuing and interest intemperature, finding newpossibly methodsbecause to im- these coagulation prove the control of conditions facilitated the overlap of these mentioned phenomena. milk coagulation. Figure 10. Figure Effect 10.Effect ofof temperature temperature andand protein protein concentration concentration on apparition on apparition of a shoulder of a shoulder in the in the curve of diffuse reflectance ratio versus time. (a) Data corresponded to skim milk adjusted to 30 g pro-g curve of diffuse reflectance ratio versus time. (a) Data corresponded to skim milk adjusted to 30 protein/kg tein/kg and and coagulated coagulated at 20, at 20, 25, 25, 30,30, 35,35, and and 4040°C◦ C andatatconstant and constantcalcium calciumchloride chlorideand and enzyme enzyme concentration. concentration.(b) (b)Data Datacorresponded corresponded to to skim milk skim adjusted milk to 30, adjusted to 50, 30, and 70 g protein/kg 50, and and coag- 70 g protein/kg and ulated at 20 °C and◦ at constant calcium chloride and enzyme concentration. Reproduced coagulated at 20 C and at constant calcium chloride and enzyme concentration. Reproduced with from [29]. permission from [Castillo, M.Z.; Payne, F.A.; Hicks, C.L.; Laencina, J.S.; López, M.-B.M.], [Modelling 4. Summary Casein Aggregation and Curd Firming in Goats’ Milk from Backscatter of Infrared Light.]; published by [J.The Dairystudy Res.], and modeling of the aggregation of enzymatically hydrolyzed casein mi- [2003]. celles is essential to have a better understanding of this stage of milk coagulation. This would The possible allow us to aggregation mechanisms achieve a better control of involved in casein the production gels are now being of enzymatically stud- coagulated ied using computational dairy products. simulations of Brownian dynamics [5,69]. With this type of approach, the role of repulsive barriers and different bonds types, as well The models formulated to date for the aggregation phase have some limitations and as the behav- ior/movements/collisions some of the [31], advantages. Holter’s model particle on theby modified aggregation mechanisms Foltmann [32], assumes andthat kinetics, there is can be studied and better understood. In other words, this technique is no overlap between proteolysis and the aggregation phase, and that the degree of prote-very useful in our understanding olysis of how is always the same particle at the time ofinteractions coagulation. influence gelation In the model properties proposed [5]. New by Payens [35– techniques for monitoring milk coagulation are also being tested, such as photon density wave spectroscopy (PDW) [70], fast Fourier transform continuous cyclic voltammetry [71], and fluorescence [72], demonstrating the continuing interest in finding new methods to improve the control of milk coagulation. Foods 2022, 11, 1243 12 of 15 4. Summary The study and modeling of the aggregation of enzymatically hydrolyzed casein mi- celles is essential to have a better understanding of this stage of milk coagulation. This would allow us to achieve a better control of the production of enzymatically coagulated dairy products. The models formulated to date for the aggregation phase have some limitations and some advantages. Holter’s model [31], modified by Foltmann [32], assumes that there is no overlap between proteolysis and the aggregation phase, and that the degree of proteolysis is always the same at the time of coagulation. In the model proposed by Payens [35–37], which was one of the first based on an aggregation kinetic, it was observed that the rate constant for the aggregation of hydrolyzed micelles is very low compared to the constants found experimentally. Dalgleish’s step function model [43], which was one of the first models to mention that a certain degree of hydrolysis of casein micelles needs to be reached before they can aggregate, is very complex, and simplifications have to be made to obtain analytical results. The step function model does not explain why a certain degree of hydrolysis of the casein micelles has to be achieved before they aggregate. In an attempt to solve this problem, Darling and van Hooydonk [45] proposed the energy barrier model, which has been extensively tested. However, it requires a simplification of Michaelis kinetics to allow calculations. The functionality model of Hyslop and Qvist [46,47], which hypothesized that the rate constant of aggregation depends on the number of active free sites in casein micelles, is based on the assumption that the aggregating particles are mostly linear, with many functional points. The adhesive hard sphere model of de Kruif [49,50], which describes the experimental data well until the end point of aggregation, does not model the hardening stage, so it is mostly related to the proteolysis reaction and the initial phase of micellar aggregation. Regarding the fractal aggregation model of Bremer et al. [55–57], it has to be assumed that the value of the fractal dimension is constant throughout the aggregation stage. However, when the gel is formed, this fractal structure tends to disappear. As a result, the model does not provide information about the next stage, the gel hardening. It has been observed in some cases that these models are formulated to model micellar aggregation up to a certain point, while others have difficulties in extrapolating their results to other, more concentrated systems. However, light scattering model proposed by Castillo et al. [29] could differentiate the aggregation from the hardening stages. 5. Conclusions Coagulation is a critical step for cheesemaking and it cannot be controlled without knowing its mechanisms and the stages that constitute it. However, despite the large number of studies that have been carried out over almost ten decades with the aim of modeling the gel assembly of enzymatically hydrolyzed casein micelles, this process is not yet known or understood in sufficient detail. On this regard, the main descriptive models for rennet-induced coagulation have been reviewed to evaluate the controversy on the existence of different phases, i.e., κ-casein hydrolysis, destabilized micellar aggregation, and gel hardening. Some authors have recently reviewed the micellar aggregation models, but without differentiating them from hardening. Indeed, most of the proposed models fail to differentiate micellar aggregation and gel hardening as two distinct but, to a certain extent, overlapping processes which are jointly responsible for gel assembly. A greater research effort is required to obtain models able to reconcile the stages of milk coagulation in a combined model, which integrates the enzymatic hydrolysis of κ-casein and the two stages constituting the non-enzymatic phase of gel formation: aggregation and hardening. Additionally, a successful model should be practical enough for industrial application so that its usage is suitable to anticipate the coagulation progress in the cheese plant for control improvement. Without a theoretical model such as this, the determination of an optimum cutting time during the manufacturing of rennet-induced cheeses will remain empiric. Foods 2022, 11, 1243 13 of 15 Author Contributions: Conceptualization, D.S. and M.C.; investigation, D.S. and Y.A.; writing— original draft preparation, D.S. and M.C.; writing—review and editing, M.C. and A.Z.; project administration, D.S. and M.C.; and funding acquisition, D.S. and M.C. All authors have read and agreed to the published version of the manuscript. Funding: This publication has been funded by V Call for Science and Technology Projects at the National University of Trujillo with Public Resources from the Canon 2021 (D.S.). This research has been funded by Generalitat de Catalunya (Project VALTEC; VALTEC13-1-0141-00) (M.C.) and by the Peruvian National Scholarship Program (RJ-145-2017-MINEDU/VMGI-PRONABEC-OBPOST) (D.S.). Institutional Review Board Statement: Not applicable. Informed Consent Statement: Not applicable. Data Availability Statement: No new data were created or analyzed in this study. Data sharing is not applicable to this article. Acknowledgments: The authors thank the Departament de Ciència Animal i dels Aliments and the Centre d´Innovació, Recerca i Transferència en Tecnologia dels Aliments (CIRTTA) for the usage of their facilities. Conflicts of Interest: The authors declare no conflict of interest. References 1. Gaygadzhiev, Z.; Massel, V.; Alexander, M.; Corredig, M. Addition of Sodium Caseinate to Skim Milk Inhibits Rennet-Induced Aggregation of Casein Micelles. Food Hydrocoll. 2012, 26, 405–411. [CrossRef] 2. Fox, P.F.; McSweeney, P.L.H. Rennets: Their role in milk coagulation and cheese ripening. In Microbiology and Biochemistry of Cheese and Fermented Milk, 2nd ed.; Law, B.A., Ed.; Chapman & Hall: London, UK, 1997; pp. 1–49. 3. Hyslop, D. Enzymatic coagulation of milk. In Advanced Dairy Chemistry-L Proteins; Fox, P., McSweeney, P., Eds.; Kluwer Academic: New York, NY, USA; Plenum Publishers: New York, NY, USA, 2003; pp. 839–878. 4. Horne, D.S.; Banks, J. Rennet-induced Coagulation of Milk. In Cheese Chemistry, Physics and Microbiology; Fox, P., McSweeney, P., Cogan, T., Guinee, T., Eds.; Elsevier Ltd.: London, UK, 2004; pp. 47–70. 5. Lucey, J. Rennet-Induced Coagulation of Milk. In Encyclopedia of Dairy Sciences; Fuquay, J., Fox, P., McSweeney, P., Eds.; Elsevier: Amsterdam, The Netherlands, 2011; pp. 579–584. 6. Chen, C.-C.; Chen, L.-Y.; Li, W.-T.; Chang, K.-L.; Kuo, M.-I.; Chen, C.-J.; Hsieh, J.-F. Influence of Chymosin on Physicochemical and Hydrolysis Characteristics of Casein Micelles and Individual Caseins. Nanomaterials 2021, 11, 2594. [CrossRef] [PubMed] 7. Castillo, M. Cutting time prediction methods in cheese making. In Encyclopedia of Agricultural, Food, and Biological Engineering; Taylor & Francis Group: Boca Raton, FL, USA, 2006; pp. 1–7. 8. Fox, P.F.; Guinee, T.P.; Cogan, T.M.; McSweeney, P.L.H. Fundamentals of Cheese Science, 2nd ed.; Springer: Boston, MA, USA, 2017; pp. 184–229. 9. Payne, F.A.; Castillo, M. Light Backscatter Sensor Applications in Milk Coagulation. In Encyclopedia of Agricultural, Food, and Biological Engineering; Heldman, D., Moraru, C., Eds.; Taylor & Francis: New York, NY, USA, 2007; pp. 1530–1534. 10. Arango, O.; Trujillo, A.J.; Castillo, M. Influence of Fat Replacement by Inulin on Rheological Properties, Kinetics of Rennet Milk Coagulation, and Syneresis of Milk Gels. J. Dairy Sci. 2013, 96, 1984–1996. [CrossRef] [PubMed] 11. Arango, O.; Trujillo, A.J.; Castillo, M. Predicting Coagulation and Syneresis Parameters of Milk Gels When Inulin Is Added as Fat Substitute Using Infrared Light Backscatter. J. Food Eng. 2015, 157, 63–69. [CrossRef] 12. Salvador, D.; Arango, O.; Castillo, M. In-Line Estimation of the Elastic Module of Milk Gels with Variation of Temperature Protein Concentration. Int. J. Food Sci. Technol. 2019, 54, 354–360. [CrossRef] 13. Harboe, M.; Broe, M.; Qvist, K. The Production, Action and Application of Rennet and Coagulants. In Technology of Cheesemaking; Law, B., Tamine, A., Eds.; Blackwell Publishing Ltd.: Oxford, UK, 2010; pp. 98–129. 14. Horne, D.; Lucey, J. Rennet-Induced Coagulation of Milk. In Cheese: Chemistry, Physics and Microbiology; McSweeney, P., Fox, P., Cotter, P., Everett, D., Eds.; Elsevier: Boston, MA, USA, 2017; pp. 115–143. 15. Sinaga, H.; Bansal, N.; Bhandari, B. Partial Renneting of Pasteurised Bovine Milk: Casein Micelle Size, Heat and Storage Stability. Food Res. Int. 2016, 84, 52–60. [CrossRef] 16. Koutina, G.; Knudsen, J.C.; Andersen, U.; Skibsted, L.H. Influence of Colloidal Calcium Phosphate Level on the Microstructure and Rheological Properties of Rennet-Induced Skim Milk Gels. LWT—Food Sci. Technol. 2015, 63, 654–659. [CrossRef] 17. Dalgleish, D.G.; Holt, C. A Geometrical Model to Describe the Initial Aggregation of Partly Renneted Casein Micelles. J. Colloid Interface Sci. 1988, 123, 80–84. [CrossRef] 18. Hill, R.D. The Effect of the Modification of Arginine Side Chains in Casein on the Coagulation of Rennin-Altered Casein. J. Dairy Res. 1970, 37, 187–192. [CrossRef] 19. Peri, C.; Pagliarini, E.; Iametti, S.; Bonomi, F. A Study of Surface Hydrophobicity of Milk Proteins during Enzymic Coagulation and Curd Hardening. J. Dairy Res. 1990, 57, 101–108. [CrossRef] Foods 2022, 11, 1243 14 of 15 20. Iametti, S.; Giangiacomo, R.; Messina, G.; Bonomi, F. Influence of Processing on the Molecular Modifications of Milk Proteins in the Course of Enzymic Coagulation. J. Dairy Res. 1993, 60, 151–159. [CrossRef] 21. Bauland, J.; Famelart, M.H.; Bouhallab, S.; Jeantet, R.; Roustel, S.; Faiveley, M.; Croguennec, T. Addition of Calcium and Magnesium Chlorides as Simple Means of Varying Bound and Precipitated Minerals in Casein Micelle: Effect on Enzymatic Coagulation. J. Dairy Sci. 2020, 103, 9923–9935. [CrossRef] [PubMed] 22. Walstra, P.; Wouters, J.; Geurts, T. Dairy Science and Technology, 2nd ed.; Taylor & Francis Group: Boca Raton, FL, USA, 2006; pp. 34–37. 23. Corredig, M.; Salvatore, E. Enzymatic Coagulation of Milk. In Advanced Dairy Chemistry; McSweeney, P.L.H., O’Mahony, J.A., Eds.; Springer: New York, NY, USA, 2016; pp. 287–307. 24. Lauzin, A.; Bérubé, A.; Britten, M.; Pouliot, Y. Effect of PH Adjustment on the Composition and Rennet-Gelation Properties of Milk Concentrates Made from Ultrafiltration and Reverse Osmosis. J. Dairy Sci. 2019, 102, 3939–3946. [CrossRef] [PubMed] 25. Berridge, N.J. An Improved Method of Observing the Clotting of Milk Containing Rennin. J. Dairy Res. 1952, 19, 328–329. [CrossRef] 26. Lazouskaya, M.; Stulova, I.; Sõrmus, A.; Scheler, O.; Tiisma, K.; Vinter, T.; Loov, R.; Tamm, M. Front-Face Fluorimeter for the Determination of Cutting Time of Cheese Curd. Foods 2021, 10, 576. [CrossRef] [PubMed] 27. Castillo, M.; Payne, F.A.; Hicks, C.L.; Lopez, M.B. Predicting Cutting and Clotting Time of Coagulating Goat’s Milk Using Difuse Reflectance: Efect of PH, Temperature and Enzyme Concentration. Int. Dairy J. 2000, 10, 12. [CrossRef] 28. Tabayehnejad, N.; Castillo, M.; Payne, F.A. Comparison of Total Milk-Clotting Activity Measurement Precision Using the Berridge Clotting Time Method and a Proposed Optical Method. J. Food Eng. 2012, 108, 549–556. [CrossRef] 29. Castillo, M.Z.; Payne, F.A.; Hicks, C.L.; Laencina, J.S.; López, M.-B.M. Modelling Casein Aggregation and Curd Firming in Goats’ Milk from Backscatter of Infrared Light. J. Dairy Res. 2003, 70, 335–348. [CrossRef] 30. McMahon, D.J.; Brown, R.J. Milk Coagulation Time: Linear Relationship with Inverse of Rennet Activity. J. Dairy Sci. 1983, 66, 341–344. [CrossRef] 31. Holter, H. Über die Labwirkung. Biochem. Zeitschr. 1932, 255, 160–168. 32. Foltmann, B. On the enzymatic and the coagulation stages of the renneting process. In Proceedings of the 15th International Dairy Congress, London, UK, 29 June–3 July 1959. 33. Foltmann, B. The biochemistry of prorennin (prochymosin) and rennin (chymosin). In Milk Proteins: Chemistry and Molecular Biology II; McKenzie, H.A., Ed.; Academic Press: New York, NY, USA, 1971; pp. 217–254. 34. Dalgleish, D.G. A New Calculation of the Kinetics of the Renneting Reaction. J. Dairy Res. 1988, 55, 521–528. [CrossRef] 35. Payens, T. On the enzyme—Triggered clotting of casein; A preliminary account. Neth. Milk Dairy J. 1976, 30, 55–59. 36. Payens, T.A.J.; Wiersma, A.K.; Brinkhuis, J. On enzymatic clotting processes I. kinetics of enzyme-triggered coagulation reactions. Biophys. Chem. 1977, 6, 253–261. [CrossRef] 37. Payens, T.A.J.; Wiersma, A.K. On enzymic clotting processes v. Rate equations for the case of armtrary rate of production of the clotting specxes. Biophys. Chem. 1980, 11, 137–146. [CrossRef] 38. McMahon, D.J.; Brown, R.J. Enzymic Coagulation of Casein Micelles: A Review. J. Dairy Sci. 1984, 67, 919–929. [CrossRef] 39. Dalgleish, D. The Enzymatic Coagulation of Milk. In Cheese: Chemistry, Physics and Microbiology; Fox, P.F., Ed.; Springer: Boston, FL, USA, 1993; pp. 69–100. 40. Dalgleish, D.G.; Brinkhuis, J.; Payens, T.A.J. The Coagulation of Differently Sized Casein Micelles by Rennet. Eur. J. Biochem. 1981, 119, 257–261. [CrossRef] [PubMed] 41. Dalgleish, D.G. Coagulation of Renneted Bovine Casein Micelles: Dependence on Temperature, Calcium Ion Concentration and Ionic Strength. J. Dairy Res. 1983, 50, 331–340. [CrossRef] 42. Brinkhuis, J.; Payens, T.A. The influence of temperature on the flocculation rate of renneted casein micelles. Biophys. Chem. 1984, 19, 75–81. [CrossRef] 43. Dalgleish, D.G. Effect of Milk Concentration on the Rennet Coagulation Time. J. Dairy Res. 1980, 47, 231–235. [CrossRef] 44. Dalgleish, D.G. A Mechanism for the Chymosin-Induced Flocculation of Casein Micelles. Biophys. Chem. 1980, 11, 147–155. [CrossRef] 45. Darling, D.F.; Van Hooydonk, A.C.M. Derivation of a Mathematical Model for the Mechanism of Casein Micelle Coagulation by Rennet. J. Dairy Res. 1981, 48, 189–200. [CrossRef] 46. Hyslop, D.B. Enzyme-Induced Coagulation of Casein Micelles: A Number of Different Kinetic Models. J. Dairy Res. 1993, 60, 517–533. [CrossRef] 47. Hyslop, D.B.; Qvist, K.B. Application of Numerical Analysis to a Number of Models for Chymosin-Induced Coagulation of Casein Micelles. J. Dairy Res. 1996, 63, 223–232. [CrossRef] 48. Stockmayer, W.H. Theory of Molecular Size Distribution and Gel Formation in Branched-Chain Polymers. J. Chem. Phys. 1943, 11, 45–55. [CrossRef] 49. De Kruif, C.G. Casein Micelles: Diffusivity as a Function of Renneting Time. Langmuir 1992, 8, 2932–2937. [CrossRef] 50. De Kruif, C.G. Casein Micelle Interactions. Int. Dairy J. 1999, 9, 183–188. [CrossRef] 51. López, M.B.; Jordán, M.J.; Granados, M.V.; Fernández, J.C.; Castillo, M.; Laencina, J. Viscosity Changes during Rennet Coagulation of Murciano-Granadina Goat Milk. Int. J. Dairy Technol. 1999, 52, 102–106. [CrossRef] Foods 2022, 11, 1243 15 of 15 52. Walstra, P.; Bloomfield, V.A.; Wei, G.J.; Jenness, R. Effect of Chymosin Action on the Hydrodynamic Diameter of Casein Micelle. Biochim. Biophys. Acta 1981, 669, 258–259. [CrossRef] 53. Horne, D.S.; Davidson, C.M. Direct Observation of Decrease in Size of Casein Micelles during the Initial Stages of Renneting of Skim Milk. Int. Dairy J. 1993, 3, 61–71. [CrossRef] 54. Horne, D.S. Casein Micelles as Hard Spheres: Limitations of the Model in Acidified Gel Formation. Colloids Surf. A Physicochem. Eng. Asp. 2003, 213, 255–263. [CrossRef] 55. Bremer, L.G.B.; van Vliet, T.; Walstra, P. Theoretical and Experimental Study of the Fractal Nature of the Structure of Casein Gels. J. Chem. Soc. Faraday Trans. 1 Phys. Chem. Condens. Phases 1989, 85, 3359. [CrossRef] 56. Bremer, L.G.B.; Bijsterbosch, B.H.; Schrijvers, R.; van Vliet, T.; Walstra, P. On the Fractal Nature of the Structure of Acid Casein Gels. Colloids Surf. 1990, 51, 159–170. [CrossRef] 57. Bremer, L.G.B.; Bijsterbosch, B.H.; Walstra, P.; van Vliet, T. Formation, Properties and Fractal Structure of Particle Gels. Adv. Colloid Interface Sci. 1993, 46, 117–128. [CrossRef] 58. Horne, D.S. Determination of the Fractal Dimension Using Turbidimetric Techniques. Application to Aggregating Protein Systems. Faraday Discuss. Chem. Soc. 1987, 83, 259. [CrossRef] 59. Worning, P.; Bauer, R.; Øgendal, L.; Lomholt, S. A Novel Approach to Turbidimetry of Dense Systems: An Investigation of the Enzymatic Gelation of Casein Micelles. J. Colloid Interface Sci. 1998, 203, 265–277. [CrossRef] 60. Lehner, D.; Worning, P.; Fritz, G.; Øgendal, L.; Bauer, R.; Glatter, O. Characterization of Enzymatically Induced Aggregation of Casein Micelles in Natural Concentration Byin SituStatic Light Scattering and Ultra Low Shear Viscosimetry. J. Colloid Interface Sci. 1999, 213, 445–456. [CrossRef] 61. Vétier, N.; Banon, S.; Ramet, J.-P.; Hardy, J. Hydratation des micelles de caséine et structure fractale des agrégats et des gels de lait. Lait 2000, 80, 237–246. [CrossRef] 62. Vétier, N.; Banon, S.; Chardot, V.; Hardy, J. Effect of Temperature and Aggregation Rate on the Fractal Dimension of Renneted Casein Aggregates. J. Dairy Sci. 2003, 86, 2504–2507. [CrossRef] 63. Walstra, P. Physical Chemistry of Foods; Marcel Dekker: New York, NY, USA, 2003; pp. 425–486. 64. Horne, D.S. Formation and Structure of Acidified Milk Gels. Int. Dairy J. 1999, 8, 261–268. [CrossRef] 65. Mellema, M.; Heesakkers, J.W.M.; van Opheusden, J.H.J.; van Vliet, T. Structure and Scaling Behavior of Aging Rennet-Induced Casein Gels Examined by Confocal Microscopy and Permeametry. Langmuir 2000, 16, 6847–6854. [CrossRef] 66. Mellema, M.; Walstra, P.; van Opheusden, J.H.J.; van Vliet, T. Effects of Structural Rearrangements on the Rheology of Rennet- Induced Casein Particle Gels. Adv. Colloid Interface Sci. 2002, 98, 25–50. [CrossRef] 67. Bauer, R.; Hansen, M.; Hansen, S.; Øgendal, L.; Lomholt, S.; Qvist; Horne, D. The structure of casein aggregates during renneting studied by indirect Fourier transformation and inverse Laplace transformation of static and dynamic light scattering data, respectively. J. Chem. Phys. 1995, 103, 2725–2737. [CrossRef] 68. Lomholt, S.B.; Worning, P.; Øgendal, L.; Qvist, K.B.; Hyslop, D.B.; Bauer, R. Kinetics of the Renneting Reaction Followed by Measurement of Turbidity as a Function of Wavelength. J. Dairy Res. 1998, 65, 545–554. [CrossRef] 69. Smykov, I.T. Kinetic of milk gelation. Part I Coagulation mechanism. In Rheology: Principles, Applications and Enviromental Impacts, 1st ed.; Karpushkin, E., Ed.; Nova Science Piblisher, Inc.: New York, NY, USA, 2015; pp. 65–82. 70. Hass, R.; Munzke, D.; Vargas Ruiz, S.; Tippmann, J.; Reich, O. Optical Monitoring of Chemical Processes in Turbid Biogenic Liquid Dispersions by Photon Density Wave Spectroscopy. Anal. Bioanal. Chem. 2015, 407, 2791–2802. [CrossRef] [PubMed] 71. Shayeh, J.S.; Sefidbakht, Y.; Siadat, S.O.R.; Niknam, K. Continuous Fast Fourier Transforms Cyclic Voltammetry as a New Approach for Investigation of Skim Milk K-Casein Proteolysis, a Comparative Study. Int. J. Biol. Macromol. 2017, 103, 972–977. [CrossRef] [PubMed] 72. Panikuttira, B.; Payne, F.A.; O’Shea, N.; Tobin, J.T.; O’Callaghan, D.J.; O’Donnell, C.P. Investigation of an in-line prototype fluores- cence and infrared backscatter sensor to monitor rennet-induced coagulation of skim milk at different protein concentrations. Int. J. Food Sci. Technol. 2020, 55, 175–182. [CrossRef]