Research and application of a new multilevel fuzzy comprehensive evaluation method for cold stress in dairy cows

Effective and comprehensive evaluation of cold stress is critical for healthy dairy cow breeding in the winter. Previous studies on dairy cow cold stress have considered thermal environmental factors but not physiological factors or air quality. Therefore, this study aimed to propose a multilevel fuzzy comprehensive evaluation (FCE) method for cold stress in dairy cows based on the analytic hierarchy process (AHP) and a genetic algorithm (GA). First, the AHP was used to construct an evaluation index system for cold stress in dairy cows from 3 dimensions: thermal environment (tempera-ture, relative humidity, wind speed, and illumination), physiological factors (respiratory rate, body surface temperature), and air quality [NH 3 , CO 2 , inhalable particulate matter (PM 10 )]. Second, the consistency test of the judgment matrix was transformed into a nonlinear constrained optimization problem and solved using the GA. Next, based on fuzzy set theory, the comment set and membership function were established to classify the degree of cold stress into 5 levels: none, mild, moderate, high, and extreme. Then, the degree of cold stress in cows was obtained using multilevel fuzzy comprehensive judgment. To investigate the effect of illumination indicators on cold stress in dairy cows, 24 prelactation cows from the south and north sides were selected for a 117-d comprehensive cold stress evaluation. The results showed that the mean mild cold stress durations were 605.3 h (25.22 d) and 725.5 h (30.23 d) and the moderate cold stress durations were 67.2 h (2.8 d) and 96 h (4.0 d) on the south and north sides, respectively. Simultaneously, generalized linear mixed model showed that there were significant correlations between the daily cold stress duration and milk yield, feeding time, lying time, and active steps in the cows on both sides. This method can reasonably indicate cow cold stress conditions and better guide cold protection practices in actual production.


INTRODUCTION
Cold stress is a dysfunctional and defensive response exhibited by animals exposed to cold environments (Brouček et al., 1991;Chase, 2020). During exposure to extreme cold temperatures, animals must generate more metabolic heat through their own body regulation to maintain their core temperature (Silva, 2012). Dairy cows are thermostatic animals and diffuse heat externally primarily through evaporation, conduction, convection, and radiation (Collier et al., 1982;Chase, 2020). When the heat production and heat dissipation of cows reach a relative balance, they are healthy, and the productive performance is normal. However, when the external environment is too cold and heat dissipation is greater than heat production, this heat balance will be disrupted (Schnier et al., 2003). This condition can lead to hormonal changes in the cow's body, causing physiological malfunctions and ultimately causing cold stress (McDowell et al., 1976;Mader et al., 2010).
Cows under cold stress have decreased production performance, reduced disease resistance, and direct disease in severe cases (Young, 1981). Wei and Lin (2010) found that severe colds had a dramatic effect on cow production and milk yield performance in North China, with milk yield decreasing by up to 40%. According to relevant findings in Harbin and Beijing, cold stress milk yield losses accounted for 8.3 and 7.5% of milk yield during the whole lactation period, respectively (Liu et al., 2009). In China, more than 68% of dairy farms are located in the northern region (Xue, 2020). The cold and long winters in these areas prevent cows from having an optimal production, and they must be given adequate attention (Zheng and Shi, 2016).
The thermal environment is an important factor affecting the cold stress of dairy cows, and the main aspects contributing to cold stress in dairy cows include Research and application of a new multilevel fuzzy comprehensive evaluation method for cold stress in dairy cows X. Fu, 1 Y. Zhang, 1 Y. G. Zhang, 2 Y. L. Yin, 1 S. C. Yan, 1 Y. Z. Zhao, 3 and W. Z. Shen 1 * temperature, humidity, wind speed, and illumination (Sun et al., 2021). Currently, researchers have proposed various cold stress evaluation methods based on the relationship between thermal environmental factors and productive performance or behavior. The temperaturehumidity index (THI) is the most widely used thermal environment index and contains both temperature and humidity parameters (Yan et al., 2019). Although this index is widely used for heat stress evaluation, it can also be used for cold stress evaluation. Xu et al. (2015) classified the THI cold stress threshold as follows: THI > 38, no stress; 25 < THI ≤ 38, mild stress; 8 < THI ≤ 25, moderate stress; −12 < THI ≤ 8, high stress; −25 < THI ≤ −12, extreme stress; THI ≤ −25, dangerous stress. Subsequently, the US and Canadian Meteorological Centers jointly proposed an updated wind chill temperature (WCT). This index is based on the temperature and wind speed and is used to assess the tolerance of humans and livestock for cold conditions (Tew et al., 2002). The specific classification criteria for the degree of cold stress are as follows: WCT > −10, no stress; −25 < WCT ≤ −10, mild stress; −45 < WCT ≤ −25, moderate stress; −59 < WCT ≤ −45, high stress; and WCT ≤ −59, extreme stress. Additionally, some studies have provided tables of the WCT of cows corresponding to different temperatures and wind speeds (Angrecka and Herbut, 2015). Dong et al. (2013) used the WCT to analyze the data from the last decade in Beijing, and the results showed that on average, 33 to 71 d per year were under mild cold stress. Temperature, humidity, and wind speed are the main drivers of heat exchange between the environment and animals. However, Webster (1971) and Keren and Olson (2006) used physical models to show that solar radiation is critical for the heat balance of cows in cold environments. Regarding lactation, 16 h of illumination per day (114-207 lx) increased the milk yield by 10 to 15% compared with a natural length photoperiod (39-93 lx) of 9 to 12 h (Peters et al., 1978). Mader et al. (2010) combined the temperature, relative humidity, wind speed and direct solar radiation time to construct a comprehensive climate index (CCI). When the CCI is less than 0, the livestock shows a cold stress response, which is classified as follows: −10 < CCI ≤0, mild stress; −20 < CCI ≤ −10, moderate stress; −30 < CCI ≤ −20, high stress; −40 < CCI ≤ −30, extreme stress; CCI ≤ −40, dangerous stress.
Currently, with the development of modern information technology, researchers have accurately obtained various physiological and behavioral characteristics of cows using wearable or contactless monitoring techniques, whereas finding that physiological responses in cows are significantly associated with cold stress (Rutten et al., 2013;He et al., 2016). Among them, the respiratory rate (RR) and ventral surface temperature (VST) are the core factors used to evaluate the internal balance of ruminant organisms and most directly reflect the health status of cows (Gonyou et al., 1979;Wang et al., 2017). Animals under cold stress can reduce the heat dissipation of the body by reducing the RR and increasing the depth of respiration. Bai and Luan (2015) found that the normal RR of dairy cows is 18 to 28 breaths per minute, which reduces to 10 to 15 breaths per minute when the temperature decreases from 10 to −15°C. The surface temperature of the skin, as the interface between the organism and environment, can reflect the changes made by the organism to adapt to the environment (Hoffmann et al., 2013). Infrared thermographs and body cooling data have shown that the ventral surface is a crucial conductive heat-loss portal in animals (Korhonen and Harri, 1986). Chen (2014) found that when high-yielding cows were in a cold breeding environment, the temperature and humidity ranged from −3.2 to 11.6°C and 70 to 96%, respectively. Their RR and VST were between 14 and 40 breaths per minute and 15.75 to 27.26°C, respectively, decreasing below normal production conditions. This finding that further analysis of physiological responses can increase the validity of the evaluation of the cold stress degree in cows.
In recent years, an increasing number of farms have chosen to keep cows in farmhouses all year round, and high-density farming has become a common form of intensive production (Winsten et al., 2010;March et al., 2014). The increase in stocking density has led to a corresponding increase in feed input and manure production per unit area, which has also considerably affected air quality in the farmhouse (Mandel et al., 2016). Manure storage, surplus feed, and individual cows produce various air pollutants, including the gaseous pollutants CO 2 , NH 3 , H 2 S, and CH 4 , and the solid pollutant particulate matter (PM; Ni, 2015). The CO 2 , NH 3 , and PM concentrations are the highest and have the largest effect on the performance and health status of dairy cows, findings that have attracted extensive attention from researchers (Cambra-López et al., 2010;Sanchis et al., 2019). Carbon dioxide is mainly derived from intestinal fermentation, respiration, and manure management, and high levels will lead to adverse effects such as chronic hypoxia and reduced appetite in cows (Jentsch et al., 2009). Ammonia is the most hazardous toxic gas to the health of dairy cows (Sanchis et al., 2019) and is mainly derived from cow manure and feed spoilage. When the concentration of NH 3 in the air rises, it can cause a high NH 3 concentration in the blood of animals, affecting the metabolism of brain nerves and muscle cells and inhibiting the feeding center (Wu et al., 2012). Particulate matter is a widespread air pollutant that comprises solid and liquid particles suspended in the air; inhalable particulate matter (PM 10 ) often carries many bacteria, viruses, and other harmful substances, such as heavy metals and volatile organic compounds. Prolonged exposure to high PM concentrations predisposes cows to respiratory and cardiovascular diseases (Cambra-López et al., 2010). When thermoregulatory mechanisms cannot compensate for cold conditions, cows initiate a series of nonspecific autoimmune responses that place the body's immune function in a suppressed state to mitigate the effects of cold stress (Marcé et al., 2010). However, various diseases caused by pollutants can lead to decreased immunity and weakening of cows (Buczinski et al., 2021), reducing their ability to withstand cold stress. Additionally, most dairy barns are built with windowed and airtight structures in the winter to improve thermal insulation, and achieving air exchange inside and outside the house only with gable and bell roof ventilation is challenging (Jungbluth et al., 2001). Inadequate ventilation in the winter often leads to higher levels of contaminants in the house that are more harmful to the health of cows and further deepens the degree of cold stress in the winter (Deng et al., 2014;Rong et al., 2014). Therefore, we also used air quality as an indicator to assess cold stress in dairy cows in this study.
From the currently proposed cold stress evaluation methods and results of related literature, the degree of cold stress in dairy cows is not a definite value but a fuzzy concept within a certain range. However, the more severe the degree of cold stress, the greater the effect on the productive performance and health status of dairy cows. Currently, the fuzzy comprehensive evaluation (FCE) method based on fuzzy mathematics is suitable to synthesize the analysis results based on the single evaluation of the target in all directions and multiple angles. In a multiparameter decision problem, each indicator must be assigned separate response weights because of their different levels of importance. The analytic hierarchy process (AHP) used in this study is a combined quantitative and qualitative indicator weight determination method with high heterogeneous data integration capability , which decomposes complex problems into several levels and several indicators by analyzing the relationships among indicators and constructing judgment matrices. The AHP is also supplemented by the expert evaluation method to assign weights to each index factor systematically and logically. The expert evaluation method is a qualitative analysis method that relies on the knowledge and experience of experts in related fields to judge, evaluate, and predict problems by ex-perts. This method has been combined with multilevel analysis or fuzzy algorithms to successfully solve specific agricultural problems. For example, Ruben et al. (2019) proposed an expert system based on a fuzzy logic algorithm whose goal is to model 6 variables (temperature, disease prevention measures for cows, annual rainfall, cow breed, production system, and breeding plan) that affect livestock production performance by applying fuzzy logic theory. These variables were identified by an expert, who assigned different fuzzy sets to each variable with knowledge and experience learned through years of observation and practice. The results show that this expert system is a useful decision support system to help farmers find the best breeding solution to maximize milk yield. José et al. (2019) used fuzzy logic and the AHP for water quality evaluation of freshwater intensive whitefish culture ponds, where the importance of key parameters (temperature, dissolved oxygen concentration, pH, and nonionized ammonia) was determined by aquaculture experts. The study findings indicate that the work provides an effective and accurate water quality evaluation program that will help water resource management and enable sustainable fish growth and reproduction.
However, the knowledge structure and cognitive level of each expert will be different, and the constructed judgment matrix sometimes cannot meet the consistency requirement. Genetic algorithm (GA) is a parallel stochastic search optimization method that simulates the genetic mechanism of nature and biological evolution (Maulik and Bandyopadhyay, 2000). In this study, the algorithm transforms the failed judgment matrix consistency test into a nonlinear constrained optimization problem. The algorithm can fine-tune the weight vector obtained from the judgment matrix to satisfy the consistency of the judgment matrix. Thus, the index weights at all levels can be determined more objectively.
Therefore, this study integrates external environmental factors affecting cow growth, as well as internal physiological factors, (1) to achieve construction of an objective, accurate, and comprehensive cold stress evaluation method for dairy cows based on the analytic hierarchy process, a genetic algorithm, and a multilevel fuzzy comprehensive evaluation model; (2) to investigate the effect of illumination factors on cold stress in cows by selecting 12 prelactating cows on each of the south and north sides of the barn for cold stress evaluation; and (3) to analyze of the significance and correlation between cold stress duration and milk yield, feeding time, lying time, and active steps of cows to verify the feasibility of the method.

Animal Management and Data Acquisition
The experiment was performed under an experimental license from Northeast Agricultural University, Harbin, China. All the experimental procedures were conducted in accordance with the principles and responsibilities outlined in the university's guidelines for animal research. The experiments were performed at a commercial farm in Harbin from November 2, 2020, to February 26, 2021. The average temperature in Harbin in the winter is −13.6°C, with extreme temperatures reaching −40°C. The lactation house is 122 m long, 30 m wide, and 5.6 m high. The experimental barn faces south, with sliding windows on the south and north sides; all windows are closed in winter. Two doors are located on the north wall, connecting the milking parlor passage in the middle, and one door is located on each of the east and west walls. The roof is a double-slope structure, using composite color steel plates and sunlight panels. A feeding channel is located in the middle of the barn. There are 287 Holstein cows in the barn. To explore the effect of illumination difference on dairy cows, we selected 24 healthy high-yielding dairy cows with 1 to 2 parity and lactation days between 24 and 80 d, with 12 cows each on north and south sides, milking 3 times a day (0500, 1100, 1900 h). The average weight of the cows was 540 ± 40 kg, and all BCS ranged from 2.75 to 3.5 (BCS is scored on a 5-point scale with each 0.25 point being an incremental unit, with 1.00 indicating a severely thin cow, and 5.00 indicating an overly fat cow). For data recording purposes, the 24 cows were specially tagged and kept in separate enclosures on the north and south sides. During the experiment, the cows were fed pellets and sheep grass twice a day (0700 and 1500 h) at a 3:7 ratio. The nutritional composition of the pellets was 49% corn, 15% bran, 11% soybean meal, 10% vegetable meal, 10% cotton meal, and 5% premix. After cutting, sheep grass accounted for 48% of the lengths greater than 19 mm, 22% of lengths 8 to 19 mm, 18% of lengths 1.2 to 8 mm, and 12% of lengths less than 1.2 mm. Figure 1 outlines the data collection points and collection process. The collection location on the south side was above the fence at a height of 1.5 m equal to that of the cattle when standing, and the north side was symmetrically arranged. The temperature, relative humidity, illumination, CO 2 , NH 3 , and PM 10 were automatically recorded using a multifunctional measuring instrument for environmental parameters (ZX-AQI-L, Shenzhen Zhongxing Environmental Protection Technology Ltd. Co.) with measurement accuracies of ± 0.1°C, ± 2%, ± 5 lx, ± 3 ppm, ± 0.5 ppm and ± 5 µg/m 3 , respectively. A thermal anemometer (GM8903, Shenzhen Jumaoyuan Technology Ltd. Co.) was used to measure the wind speed around the cattle body 3 times continuously, and the resolution was 0.01 m/s. A handheld laser infrared thermometer (HT-8962, Shanghai Hongcheng Technology Ltd. Co.) was used to determine the VST of the cows. The measurement was continued 3 times, with a measurement accuracy of ± 0.1°C. The farm staff obtained the RR by observing the rise and fall of the cow's venter, and the value was recorded 3 times continuously for 1 min each time. A multifunctional measuring instrument for environmental parameters can set the collection time interval. The thermal anemometers, laser infrared thermometers, and manual measurements of cow RR were performed by 4 trained researchers and 2 farm managers. Before the start of the experiment, 6 people practiced monitoring the wind speed, RR, and VST, together for 2 d to ensure the similarity of each data monitoring procedure. In this study, we needed to obtain daily milk yield, so we defined the data collection cycle as lasting from the start of the previous day last milking (1900 h) to the last milking of the next day (1900 h). Data were collected simultaneously on the south and north sides before milking with a 2-h collection interval, and a total of 33,696 data points were obtained. This period was also used as the period for calculating the daily cold stress duration in cows.
A total of 117 d of experimental data collection were performed on 24 cows. The temperature ranged from −9.9 to 6.9°C and −10.6 to 6.6°C, and the humidity ranged from 34.8 to 100% and 41.3 to 99.5%, on the north and south sides, respectively. The wind speed ranged from 0 to 0.32 m/s on both sides. The illumination ranged from 0 to 869 lx and 0 to 442 lx, and the illumination times were 17.1 and 16 h/d, respectively; artificial illumination intensity was 79 lx. The CO 2 concentration ranged from 367 to 4,401 ppm and 508 to 4,965 ppm, the NH 3 concentration ranged from 0 to 8.9 ppm and 0 to 7.9 ppm, and the PM 10 concentration ranged from 23 to 468 µg/m 3 and 24 to 463 µg/ m 3 on the north and south sides, respectively. The RR and VST were the average values of 12 cows on the north and south sides at each collection moment. The RR ranged from 16.28 to 29.43 breaths per minute and 16.30 to 29.1 breaths per minute. The VST ranged from 12.16 to 33.97°C and 12.03 to 33.90°C on the north and south sides, respectively. Figure 2 shows the daily averages of the 9 index factors.
Additionally, to investigate the relationship between cold stress duration and the behavioral characteristics of cows, 3 types of behavioral data (including feeding time, lying time, and active steps) were collected using a smart collar and smart leg ring developed by RuiBaoLe Dairy Division. The device works by monitoring the change in the angle between the smart collar and smart leg ring to identify the behavioral characteristics of the cow. The device transmits data to the server every 5 min, and the server records the daily behavior data by accumulating them. The daily milk yield data were collected and provided by the farm. To ensure the accuracy of the data, the wearable device was installed 7 d before the experiment.

Data Processing and Analysis
Because wind speed, RR, and VST were measured at 3 time points in this study, we can determine the measured value of cold stress at these 3 different times via the cold stress evaluation model. These 3 measurements were correlated, and RR, VST, milk production, and behavioral characteristics were dependent variables with non-normal distributions. Therefore, this study proposed the generalized linear mixed model (GLMM) to analyze the relationship between average daily RR, VST, and daily milk yield, and the relationship between daily cold stress duration, milk production and cow behavioral characteristics. The GLMM is an extended form of generalized linear model (GLM) and linear mixed model (LMM; Berry and Evans, 2022) that is thought to provide unique advantages when analyzing repeated measures data (Gibbons et al., 2010). And the maximum likelihood method was used to estimate the parameters of the GLMM. The model enables greater rationality and reliability of the results, thereby confirming the validity of the comprehensive model for evaluating cold stress in dairy cattle. The GLMM model is as follows: where Y is the n × 1-dimensional vector, representing the observed values (i.e., the dependent variable); μ is the observed prediction vector; ε is a random error term; g −1 (·) is the inverse function of the monotone differentiable link function g(·); X is the N × p matrix, representing the independent variables of fixed effects; Z is the N × p matrix, representing the independent variables of random effects; and β and γ are vectors of p × 1, representing the fixed effect and the random effect parameter vectors of the model, respectively. The best-fitting model was selected based on the Akaike information criterion (AIC) and the Bayesian information criterion (BIC). The test level was set at α = The temperature, relative humidity, illumination, CO 2 , NH 3 , and inhalable particulate matter (PM 10 ) were automatically recorded using a multifunctional measuring instrument for environmental parameters. (b) A thermal anemometer was used to measure the wind speed around the cattle body. (c) A handheld laser infrared thermometer was used to determine the ventral surface temperature of the cows. (d) The smart collar and smart leg ring were used to obtain real-time data on the cow's feeding time, lying time, and activity steps. 0.05, and the P-value represents the significant level: P < 0.01 is considered highly significant, and P < 0.05 is considered significant. P-values less than 0.001 are reported as P < 0.001. The data were presented as means ± standard deviation.

The Proposed Methodology
The analytic hierarchy process, genetic algorithm, fuzzy comprehensive evaluation (AHP-GA-FCE) method proposed in this study comprises the following 4 steps ( Figure 3).
Step 1: According to the literature, the factors related to the cold stress of dairy cows are classified into characteristics, and the comprehensive evaluation index system is constructed using the AHP.
Step 2: The judgment matrix of each index layer is established, and the consistency test is performed. If the consistency requirement is satisfied, the weights of each index layer are determined. By contrast, the weights are optimized using the GA.
Step 3: The GA is mainly used to modify the judgment matrix and calculate the weights of each index layer through crossover and variation operations, and finally obtain the optimal feasible solution.
Step 4: On this basis, the comment set and membership function were set from the perspectives of the degree of cold stress and applicability, and the multilevel index factors were evaluated comprehensively. Finally, the degree of cold stress of dairy cows was determined according to the principle of maximum membership.

Construction of the Evaluation Index System
The evaluation index system was divided into 4 layers, including the target layer, constraint layer, criterion layer, and scheme layer ( Figure 4). The target layer is the core of this study. The constraint layer is the firstlevel index, which includes 3 aspects of the thermal environment, physiological factors, and air quality. The criterion layer is the 2-level index that further reflects the content of the constraint layer, and includes 9 index factors (temperature, relative humidity, wind speed, illumination, RR, VST, CO 2 , NH 3 , and PM 10 ). Finally, 5 cold stress evaluation degrees were divided into no, mild, moderate, high, and extreme cold stress in the scheme layer.

Construction of the Judgment Matrix
Because the importance of each index affecting cold stress in dairy cows varies, each index should be assigned a corresponding weight. Assuming that n factors affect the evaluation object, U = {u 1 , u 2 , …, u n } is the set of factors, and u i (i = 1, 2, ···,n) is the influencing factor. Each factor is grouped and set as U = {U 1 , U 2 , ···, U n }, where m influencing factors are contained in U i , factor, and u im is the 2-level index factor. Additionally, the following 2 conditions must be met: The relative importance of different factors in the set of factors U and U i in each index layer was quantified by experienced and knowledgeable experts or breeders using a scale of 1 to 9 (  Step 1: The factors related to the cold stress of dairy cows are classified into characteristics, and the comprehensive evaluation index system is constructed using the AHP. Step 2: The judgment matrix of each index layer is established, and the consistency test is performed. If the consistency requirement is satisfied, the weights of each index layer are determined. By contrast, the weights are optimized using the genetic algorithm (GA).
Step 3: The GA is mainly used to modify the judgment matrix and calculate the weights of each index layer through crossover and variation operations, and finally obtain the optimal feasible solution.
Step 4: The comment set and membership function were set from the perspectives of the degree of cold stress and applicability, and the multilevel index factors were evaluated comprehensively.
where a ij is the importance of factors u i and u j , and a ij > 0, a ii = 1, a ji = 1/a ij .

Calculation of the Weight Vector
In this study, the square root method was used to determine the weights of each factor. The product of each row of the judgment matrix (M i ) is calculated as follows: , , , .
The weight vector w i ( ) of each factor set was calculated as follows: The vector was normalized as follows: , , , .

Testing Consistency of the Judgment Matrix
To ensure the rationality of the evaluation results, the judgment matrix was tested for consistency. The consistency index (CI) of the judgment matrix was calculated as follows: where n is the rank of the judgment matrix, and n > 1; λ max is the maximum characteristic root of the judgment matrix; and (Aw) i denotes the ith component of Aw.
The consistency ratio (CR) of the judgment matrix was calculated as follows: where RI is the random consistency index of the judgment matrix (Xiang et al., 2014). If 0 ≤ CR < 0.1, the judgment matrix satisfies the consistency requirement and the vector w is used as the solution of the weight vector. Otherwise, if CR ≥ 0.1, the weights of the judg- . Comprehensive evaluation index system for cold stress in dairy cows. The target layer is the core of this study; the constraint layer is the first-level index; the criterion layer is the 2-level index that further reflects the content of the constraint layer, and includes 9 index factors; and, 5 cold stress evaluation degrees were divided into no, mild, moderate, high, and extreme cold stress in the scheme layer. ment matrix must be optimized to meet the required consistency.

Optimizing the Weight and Consistency of the Judgment Matrix
In this study, GA was used to optimize the consistency test problem in the judgment matrix to make the judgment matrix meet the consistency. The model was optimized using the following steps.
Step 1: Considering the judgment matrix of the firstlevel evaluation index (U = {U 1 , U 2 , ···, U n }) as an example, and its weight vector is (W = w 1 , w 2 , ···, w i ). The optimization objective function is constructed as follows: where CIF (n) is the consistency index function and w k is the optimized weight vector. The weight value corresponding to the minimum value of CIF (n) is the optimal weight value corresponding to the judgment matrix U. If CIF (n) = 0, the judgment matrix U has complete consistency.
Step 2: Random initialization of populations: The initial population is randomly set to generate a set of feasible solutions. The individual coding method selects real number coding, and each individual is a real number string. The population size is n, the number of evolutions is 200, the crossover probability is 0.3, and the mutation probability is 0.2.
Step 3: Construction of the fitness function: The sum of the absolute values of the errors of the predicted and expected outputs of CIF (n) is used as the individual fitness value F i as follows: where n is the number of network output nodes; Y i is the expected output of the ith node of the CIF (n); X i is the predicted output of the ith node; and k is the coefficient.
Step 4: Selection operation: The roulette method is chosen as the strategy for the adaptation ratio, and the selection probability of each individual i is p i .
where F i is the fitness value of individual i. Because smaller fitness values are better, the reciprocal of the fitness value is obtained before individual selection; n is the number of individuals in the population.
Step 5: Crossover operation: Using the real number crossover method, the kth chromosome (c k ) and lth chromosome (c l ) at position j are crossed over as follows: . The 1 to 9 scaling method; u i and u j are the factors affecting the evaluation object. 1 represents 2 factors of equal importance; 3 represents u i /u j is weakly important compared with u j /u i ; 5 represents u i /u j is important compared with u j /u i ; 7 represents u i /u j is extremely important compared with u j /u i ; 9 represents u i /u j is absolutely important compared with u j /u i . Values 2, 4, 6, and 8 are the intermediate values of the above judgments.
where b 1 is the crossover probability.
Step 6: Mutation operation: the jth gene (c ij ) of the ith individual is selected for mutation as follows: where c max is the last generation of c ij and c min is the next generation of c ij ; r 1 is the mutation probability; r 2 is the random number between [0,1], and r 2 = 0.2; g is the number of current iterations, and G max is the maximum number of evolutions.
Step 7: Whether the updated weight meets the constraints of step 1 is determined. If so, the iteration is terminated, and the corresponding optimal feasible solution is the optimal weight value. Otherwise, the process returns to step 3.

Construction of the Multilevel FCE Model
Setting the Comment Set. Because the evaluation value of each index factor is different, different degrees are often formed, and the set comprising different judgments is called the comment set. Based on the evaluation index system for cold stress in dairy cows ( Figure  4), the evaluation index factor set and comment set were set as follows: 1) Set the first-level evaluation index factor set U = {(U 1 , U 2 , U 3 )}, and the 2-level evaluation index factor set U 1 = (u 11 , u 12 , u 13 , u 14 ), U 2 = (u 21 , u 22 ), U 3 = (u 31 , u 32 , u 33 ). 2) Set the comment set V j = (V 1 , V 2 , ···, V j ), (j = 1, 2, ···, 5); V 1 -V 5 represent 5 cold stress degrees: no, mild, moderate, high, and extreme, respectively.

Construction of the membership of evaluation indices.
The classification of 9 index factors of cold stress degree was completed based on environmental quality standards and current literature studies, and quantitative and qualitative analyses were performed.
To ensure the quality of livestock products, largescale dairy farms in North China should follow the specified environmental quality standards of Beijing Agriculture Bureau (2007) in the actual breeding process. This standard is used for environmental quality control, monitoring and environmental management of dairy farms. The ranges and thresholds of temperature, relative humidity, wind speed, illumination, CO 2 , NH 3 , and PM 10 are shown in Table 1.
Temperature is the primary factor affecting cold stress in dairy cows. Brouček et al. (1991) and Angrecka and Herbut (2015) both found that cows perform best when they are in a breeding environment of 4 to 20°C. MacDonald and Bell (1958) and Bai and Luan (2015) showed that the milk yield of cows decreased by 5 and 10% when the temperatures were below −4°C and −20°C, respectively. Zhang (2015) also reported that the milk yield of cows decreases by more than 6% when the temperature is below −10°C. Therefore, the temperature was divided into 5 parts: 4 to 20°C, −4 to −10°C, −15 to −20°C, −25 to −30°C, and ≤−35°C, indicating no, mild, moderate, high, and extreme cold stress degrees, respectively. Ma (2017) showed that excessive humidity in the house in the winter promoted the reproduction of microorganisms and improves the probability of dairy cows developing skin diseases such as eczema, and that humidity should be controlled between 50 and 70%. Zhao et al. (2018) and Tang et al. (2019) showed that if the humidity in the house was higher than 80% in winter, it was unsuitable for the normal production of dairy cows. Jing and Jing (2021) suggest that the upper limit of humidity in dairy barns should be 85%, and when the air is nearly saturated with humidity (100%), the moist skin surface makes the cow colder. Therefore, the humidity is classified as 50 to 70%, 75 to 80%, 83 to 87%, 90 to 95%, and ≥97%.
The wind speed requirements in the cattle barn are different in the summer and winter. Qi (2013), Yu et al. (2020), and Jing and Jing (2021) believed that the  suitable wind speed range in winter was 0.1 to 0.2 m/s, not exceeding 0.25 m/s. A wind speed higher than 0.2 m/s will significantly increase the nonevaporative heat dissipation, increase the heat production of the cows, and make the cows feel cold. Below 0.1 m/s, ventilation is insufficient, which prevents the discharge of moisture and harmful gases. The wind speed was classified as 0.1 to 0.2 m/s, 0.25 to 0.3 m/s, 0.35 to 0.4 m/s, 0.45 to 0.5 m/s, and ≥0.55 m/s. The degree of cold stress must be classified in various ways for illumination factors.
In actual production, the milk yield of cows receiving supplementary illumination (50.2 kg/d) was higher than that of cows under normal illumination conditions (46.6 kg/d; Vanbaale et al., 2005). Farming conditions with artificial illumination intensity greater than 150 lx and uniform brightness in all areas of the house are more conducive to the healthy breeding of cows (Peters et al., 1978;Pajumägi et al., 2007). According to the 3 aspects of illumination duration, artificial illumination intensity and illumination distribution uniformity, experts combined their experience to qualitatively analyze the degree of cold stress of illumination indicators.
In the environmental quality standard, the CO 2 concentration threshold is 835 ppm. However, during the studies on dairy farming, Qi (2013) found that the normal breathing and health of cows are only affected when the CO 2 concentration in the house is higher than 1,500 ppm. Vtoryi et al. (2017) and Karandušovská et al. (2016) both agreed that cows can tolerate a CO 2 concentration of 2,500 ppm in the barn. Therefore, the CO 2 concentration was classified as 0 to 800 ppm, 900 to 1,500 ppm, 1,800 to 2,500 ppm, 2,800 to 3,500 ppm, and >4,000 ppm. Qi (2013) concluded that when NH 3 is in the range of 15 to 20 ppm, it reduces the resistance of animals to diseases. However, Karandušovská et al. (2016) considers the limit value of the NH 3 concentration in the dairy barn to be 25 ppm, whereas the environmental quality standards indicate that the NH 3 concentration must be less than 28 ppm. Sanchis et al. (2019) and Cambra-López et al. (2010) found that NH 3 and PM 10 are toxic gases and that their concentration must be as low as possible. Therefore, NH 3 was divided into 5 parts: 0 to 10 ppm, 13 to 15 ppm, 17 to 20 ppm, 23 to 25 ppm, >27 ppm. The PM 10 concentration in the dairy barn should not exceed 3 mg/m 3 ; PM 10 was classified as 0 to 1 mg/m 3 , 1.3 to 1.5 mg/m 3 , 1.7 to 2.0 mg/m 3 , 2.3 to 2.5 mg/m 3 , and >3 mg/m 3 .
Most indices for stress evaluation in dairy cows are constructed based on lactation (Sun et al., 2021). The milk yield of dairy cows is the most sensitive to climate change; yet, milk yield is the most crucial factor. Based on early THI studies using artificial climate chambers to investigate the effects of different temperatures and humidity on lactation in Holstein cows, the relationship between the THI and decreased lactation in cows was evaluated, and equations to estimate decreased lactation using the THI were derived (Berry et al., 1964). Linvill and Pardue (1992) calculated THI values at 1-h intervals and found that the highest correlation between total hours and daily milk production when THI values exceeded 74 or 80. And they developed a model linking total hours with THI values exceeding 74 and 80 to daily milk production. The black globe humidity index (BGHI) was based on lactation, completing a further adaptation of the index calculation model (Buffington et al., 1981). The equivalent temperature index (ETI) was also established based on milk yield (Silva, 2012). Therefore, in this study, the average daily milk yield of cows of the same lactation age in the spring was used as a benchmark, and RR and VST were classified according to the degree of milk yield reduction.
Quantitative Indices. The basic idea of the fuzzy mathematical model is the idea of membership. The membership degree of each evaluation level of cold stress in cows can be transformed by the corresponding fuzzy mathematical language. This study selected trapezoidal and semi-trapezoidal distributions to construct the membership function. Each distribution is divided into minor membership, intermediate membership, and large membership (Ruben et al., 2019), and the range of membership is 0 to 1. As shown in Figure 6, the 8 index factors (temperature, relative humidity, wind speed, RR, VST, CO 2 , NH 3 , and PM 10 ) were quantified using the corresponding membership functions with the following formulas: where μ M (x i ), μ I (x i ), and μ L (x i ) are minor, intermediate, and large membership functions respectively. Qualitative Indices. The data of this study showed that the average daily illumination intensity ranged from 0 to 869 lx and 0 to 442 lx on the north and south sides, with illumination durations of 17.1 and 16 h/d, respectively. The artificial illumination intensity was 79 lx. The uniformity of distribution is judged by how the lights are hung in the barn, distance, range of illumination, and presence of alternating light and dark in the barn. The collected data and barn monitoring video were sent to the experts to score the illumination index, and the average of experts was taken as the final score of the illumination index. Table 2 and Table 3 show the evaluation content and scoring criteria. The index score corresponds to the degree of cold stress in cows, and higher scores indicate a more positive effect of the illumination index on milk yield.
Comprehensive Evaluation. (1) Two-level fuzzy evaluation. First, the u im th of the 2-level index factor set (U i = {u i1 , u i2 , ···u im }) was judged, and the singlefactor judgment matrix was obtained. Thus, the membership of u im to V j is determined, and the results of the single-factor evaluation result are obtained: The single-factor evaluation is a fuzzy mapping of the set of index factors to the set of comments. Judging each index factor in U i , the fuzzy judging matrix R i of U i can be obtained as follows: The weight of the 2-level index factor (U i ) is W i = (w i1 , w i2 , ···, w im ), and the results of the comprehensive evaluation were obtained as follows: (2) First-level fuzzy evaluation. The 2-level fuzzy judgment set B i = (B i1 , B i2 , ···, B ij ) T is used as the fuzzy relationship matrix for the first-level comprehensive judgment.
The weight of the first-level index factor is W = (w 1 , w 2 , ···, w i ), and the vector of the comprehensive evaluation was obtained as follows: where is the synthetic operator of the fuzzy matrix; b j is the membership degree of the evaluation object for each comment; and finally, the degree of cold stress of dairy cows is obtained according to the principle of the maximum membership degree.

RESULTS
In the study, the average daily RR and VST of cows were considered dependent variables in the GLMM. The fixed effects component consisted of an intercept and the daily milk yield of cows, and the random effects component included individual cows. The fixed effects term was statistically significant (F = 1270.122, P < 0.001; F = 2319.157, P < 0.001), with intercepts of −25.54 and −92.20 (P < 0.001). The main effect of daily milk yield was significant (all P < 0.001), with coefficient values of 1.43 and 3.41 (Figure 7). For every 1 kg decrease in milk yield, the RR decreased by 1.43 breaths per minute, and VST decreased by 3.41°C. Additionally, the covariance parameter of the random effects term was not significant (P = 0.074 > 0.05 for the south side; P = 0.096 > 0.05 for the north side), indicating no significant differences between individual cows in terms of average daily RR and VST. Through consultation with experts, both RR and VST were divided into 5 degrees according to a 5% milk yield reduction. The 5 membership function charts for each of the 8 index factors are provided in Figure 8.

Index Weights of Dairy Cow Cold Stress
We conducted a questionnaire survey of 30 experts (12 experts in dairy nutrition and feed science, 13 experts in dairy production, and 5 experts in intelligent dairy farming) in the field of dairy research and 13 dairy farmers, using the 1 to 9 scale method. The importance of each level of index factors was scored, and 43 sets of weight vectors were calculated and averaged. Table 4 shows only the process of calculating the comparatively important scores of indicators at each level and their weight vectors by one expert and one farmer.
The weight values of each index factor were obtained (Figure 9).

Evaluation of the Cold Stress Duration of Dairy Cows
Consider the following data of any cow on January 10 as an example, where temperature: −7°C, relative humidity: 84%, wind speed: 0.23 m/s, illumination: 76 score, RR: 21.3 breaths per minute, VST: 24.2°C, CO 2 : 2,639 ppm, NH 3 : 5.3 ppm, PM 10 : 230 µg/m 3 . The weights of their corresponding 2-level index factors (u 11 , u 12 , u 13 , u 14 , u 21 , u 22 , u 31 , u 32 , u 33   ate and high cold stress, no cold stress, and no cold stress degree, respectively. The fuzzy evaluation matrix R 1 , R 2 , R 3 of the first-level index U 1 ,U 2 ,U 3 is obtained.   After normalization, B = [ ] 0 10 0 65 0 23 0 03 0 . , . , . , . , . The maximum b 2 was 0.65, and according to the principle of maximum membership, the comprehensive evaluation of cold stress belonged to the V 2 in the V. Therefore, the comprehensive evaluation result was "mild cold stress." The collected data show that the temperature, relative humidity, wind speed and pollutant concentration in dairy farming environment are not obviously changed by the position of the north and south sides. However, the illumination index was obviously different on both sides. The illumination intensity and average daily illumination duration on the south side were higher than those on the north side by 427 lx and 1.1 h, respectively. Because of artificial supplemental lighting (79 lx), the average daily illumination duration on both sides was not very different. But the illumination intensity on the south side was obviously higher than that on the north side. An uneven illumination distribution was observed because of the long hanging distance and small irradiation range between the fill lights. The experts scored the illumination index using an average value of 76 and 58, indicating that the cows were in mild and moderate cold stress on both sides. The results also have implications for evaluation of the degree of cold stress on both sides. Figure 10 shows the average daily cold stress duration on the south and north sides; mild cold stress was primarily observed, but of a long duration. The average total durations of mild cold stress were  and 96 h (4.0 d) on the 2 sides, observed primarily in December, January, and February. The duration of cold stress on the south side was 120.2 and 28.8 h less than that on the north side. Figure 11 shows the total cold stress duration, and average daily milk yield, feeding time, lying time, and active steps of 24 lactating cows on the north and south sides. Overall, the average daily milk yield and behavioral characteristics were affected by cold stress. As the duration of cold stress increased, the milk yield continued to decline. The milk yield on both sides were 34.4 ± 0.75 kg/d and 33 ± 0.26 kg/d, which decreased by 4.97 and 8.84%, respectively, compared with the average milk yield (36.2 ± 0.34 kg/d) of dairy cows in the same lactation period in the spring. The feeding time of cows on both sides were 6.45 ± 0.45 h/d and 7.29 ± 0.36 h/d, the lying time was 11.2 ± 0.37 h/d and 12.78 ± 0.36 h/d, and the active steps were 4,897 ± 345 times/d and 3,037 ± 513 times/d, respectively. Compared with the south side, lactating cows on the north side spent a longer time feeding and lying down and took fewer active steps.

Link Between Cold Stress Duration and Daily Milk Yield.
In this study, the daily milk yield and behavioral characteristics of cows were considered dependent variables in the GLMM. The fixed effects component consisted of an intercept and the daily cold stress duration of cows, and the random effects component included individual cows.
The productivity of dairy cows is mainly reflected in milk yield. Cold stress increases the resting metabolic rate of dairy cows, resulting in an increase in the energy required to maintain individuals and thus decreasing the efficiency of feed conversion into milk yield. The fixed effects term was statistically significant (F = 6911.687, P < 0.001 for the south side; F = 5787.778, P < 0.001 for the north side), with intercepts of 35.805 and 34.003 (all P < 0.001, Table 5), when the daily milk production was used as the dependent variable on the south and north sides. The main effect of daily cold stress duration was significant (all P < 0.001), with coefficient values of −0.234 and −0.136, respectively, indicating that milk yield decreases with increasing cold stress duration. With a 1 h increase in cold stress duration, milk production decreased by 0.234 and 0.136 kg. Although the decline in milk yield was smaller on the north side, daily milk yield was obviously lower than that on the south side. Additionally, the covariance parameter of the random effects term was not significant (P = 0.245 > 0.05 for the south side; P = 0.243 > 0.05 for the north side), indicating no significant differences between individual cows in terms of daily milk yield. Figure 12 illustrates the linear relationship between the daily cold stress duration and the daily milk yield of cows on both sides.

Link Between Cold Stress Duration and Behavioral Characteristics.
Dairy cow behavior is the external manifestation of spontaneous regulatory changes or physiological needs within the organism. Feed intake is often stimulated when cows are under cold stress because they can take in more energy by consuming more food, which can be used to produce heat to maintain body temperature changes. The fixed effects term was statistically significant (F = 3789.682, P < 0.001 for the south side; F = 5222.622, P < 0.001 for the north side), with intercepts of 6.054 and 6.830 (all P < 0.001, Table 5), when the daily feeding time was used as the dependent variable on the south and north sides. The main effect of daily cold stress duration was significant (all P < 0.001), with coefficient values of 0.070 and 0.065, respectively, indicating that feeding time increases with increasing cold stress duration. With a 1 h increase in cold stress duration, feeding time increased by 0.070 and 0.065 h. Additionally, the covariance parameter of the random effects term was not significant (P = 0.240 > 0.05 for the south side; P = 0.241 > 0.05 for the north side), indicating no significant differences between individual cows in terms of daily feeding time.
The fixed effects term was statistically significant (F = 2879.380, P < 0.001 for the south side; F = 3233.210, P < 0.001 for the north side), with intercepts of 10.462 and 12.182 (all P < 0.001, Table 5), when the daily lying time was used as the dependent variable on the north and south sides. The main effect of daily cold stress duration was significant (all P < 0.001), with coefficient values of 0.124 and 0.085, respectively, indicating that lying time increases with increasing cold stress duration. With a 1 h increase in cold stress duration, lying time increased by 0.124 and 0.085 h. Additionally, the covariance parameter of the random effects term was not significant (P = 0.240 > 0.05 for the south side; P = 0.241 > 0.05 for the north side), indicating no significant differences between individual cows in terms of daily lying time.
The fixed effects term was statistically significant (F = 2608.721, P < 0.001 for the south side; F = 2885.138, P < 0.001 for the north side), with intercepts of 5764.055 and 3540.489 (all P < 0.001, Table 5), when the daily active steps was used as the dependent variable on the north and south sides. The main effect of daily cold stress duration was significant (all P < 0.001), with coefficient values of −150.839 and −71.747, respectively, indicating that active steps decrease with increasing cold stress duration. With a 1 h increase in cold stress duration, active steps decreased by 150.839 and 71.747 units. Additionally, the covariance parameter of the random effects term was not significant (P = 0.241 > 0.05 for the south side; P = 0.240 > 0.05   for the north side), indicating no significant differences between individual cows in terms of daily active steps. Figure 13 illustrates the linear relationship between the daily cold stress duration and the daily feeding time, lying time and active steps in the cows on the north and south sides.

A New Method to Evaluate of Cold Stress in Dairy Cows
Generally, dairy cows are more tolerant of cold conditions than other livestock, thus, fewer studies have investigated cold stress in dairy cows. However, the effects of cold stress on dairy cows cannot be ignored and warrant further study. In this work, we not only used environmental indicators, including temperature, relative humidity, wind speed, illumination, NH 3 , CO 2 , and PM 10 in the dairy barn, but also integrated the physiological characteristics of the RR and VST of the cows. An AHP-GA-FCE based cold stress evaluation method for dairy cows is proposed. This method constructing a 2-by-2 comparison judgment matrix is the basis for the AHP and is performed by experienced and knowledgeable dairy experts. Simultaneously, a judgment matrix optimization model based on the GA is constructed that fully solves the inconsistency problem in the judgment matrix. Additionally, membership theory in the fuzzy comprehensive evaluation method The link between daily cold stress duration and daily milk yield. Panel A shows the link between the daily cold stress duration and daily milk yield for cows on the south side (P < 0.001); B shows the link between the daily cold stress duration and daily milk yield for cows on the north side (P < 0.001). Figure 13. The link between daily cold stress duration and daily behavioral characteristics of dairy cows. Panels A 1 -A 3 show the link between daily cold stress duration and daily feeding time, lying time, and active steps for 12 cows on the south side (P < 0.001); panels B 1 -B 3 show the link between daily cold stress duration and daily feeding time, lying time, and active steps for 12 cows on the north side (P < 0.001). is a reliable approach to transform human experience and knowledge into mathematical models (de Mol and Woldt, 2001;Ruben et al., 2019). Based on the above analysis, the present study transferred the experience of literature research and the environmental quality standards of northern winter into a fuzzy evaluation model. This comprehensive model overcomes the subjectivity of expert experience, which is challenging to quantify, and avoids ignoring index weights by fuzzy evaluation methods. This model is better than the simple index model, which relies solely on the fitting of thermal environmental factors and can reasonably reflect the cold stress condition of cows.

Effects of Environmental Factors on Cold Stress in Dairy Cows
The thermal environment is the physical environment that directly affects the heat dissipation process of the animal body and is closely related to the metabolic heat and thermoregulation process (Keren and Olson, 2006). The temperature determines the convective heat transfer temperature difference between the surface of the animal body and environment, affecting the convective heat transfer. Humidity is a parameter that indicates the water vapor content of the air, and higher humidity will increase the feeling of coldness in cattle. The wind speed also affects the diffusion coefficient of water vapor on the surface of the animal body, and a higher airflow speed will increase heat dissipation. Solar radiation not only raises the temperature in the house, but the solar radiation absorbed by livestock may be several times the amount of heat produced by their metabolism (Yan et al., 2019). A suitable thermal environment promotes good health and high lactation performance of dairy cows. In our study, the lowest temperature in the barn was −10.6°C, resulting in a decrease in milk yield of 8.84%. The reason is that the temperature in Harbin decreases suddenly in the autumn and remains low in the winter, more seriously affecting milk yield performance when cows are under long-term cold stress. The present study found that the high relative humidity in the barn, up to 100%, because of the cows' respiration, evaporation of feces and urine, and inadequate ventilation, caused stress response in the cows.
In our experiment, all the windows of the barn were closed in the winter for heat preservation. The windows were only opened for ventilation at noon and during milking, and the maximum wind speed was 0.32 m/s, which had little effect on the cold stress of cows. Li et al. (2009) showed that the illumination duration can be fully used to increase the house temperature in the winter, leading to all increased milk yields. Additionally, supplemental illumination promotes lactation, reproduction, growth, and development, and improves immunity in cows (Peters et al., 1981). Dahl et al. (2000) also showed that increasing illumination from less than 12 h/d (short daylight cycle) to 16 to 18 h/d (long daylight cycle) increased the average daily milk yield by 2.5 kg per head of cow. In the study, the illumination data from the north and south sides were collected and analyzed, and the experts scored the illumination index using average values of 76 and 58, respectively. The scores also affected the final evaluation of cold stress in cows. The average total cold stress duration of cows on the south side was 672.5 h, which was 149 h less than that on the north side (821.5 h). The cows on the south side had better illumination conditions, and their milk yield was 1.4 kg higher than that on the north side, a finding that was consistent with that reported by Peters et al. (1981). Therefore, in the actual breeding process, the cattle farm should increase the position and number of roof sun panels according to the orientation of the barn, sun altitude and azimuth, and regulate the use of fill lights. Thus, both sides of the herd could receive the maximum amount of solar radiation, resulting in an even distribution of illumination and a reduced duration of stress for the cows.
Thermal environmental parameters are the direct factors that cause cold stress in dairy cows, but good health can also enhance the ability to resist wind chills (Tucker et al., 2007). The weights of CO 2 , NH 3 , and PM 10 in the air quality indices were 0.05, 0.01, and 0.02, respectively (Figure 9). Although these weights were small, the effects of harmful gases on cows cannot be ignored. When the concentration of harmful gases and particulate matter in the barn accumulates to a certain level, it will not only affect the normal breathing and physiological metabolism of cows but also threaten the health of farm workers and cause human and animal respiratory diseases. In this work, the NH 3 and PM 10 concentrations in the lactation house did not reach cold stress levels, likely because of the hay laid on the floor of the house. Additionally, farm workers clean up manure, feed, and other waste from the barn twice a day in the morning and evening, a task that helps to reduce the concentration of NH 3 and PM 10 . However, the concentrations of CO 2 can reach 4,965 ppm in the study, which far exceeds the maximum tolerance range of 2,500 ppm for dairy cows (Karandušovská et al., 2016). Jentsch et al. (2009) also showed that the CO 2 concentration reflects the ventilation conditions and air pollutants in the barn, serving as an indirect indicator to evaluate the air conditioning, and a CO 2 concentration that is too high will affect the health of cows.
Furthermore, the degree of cold stress on environmental factors was classified based on the environ- mental quality standards of northern farming with existing studies. The data collected showed that the temperature, relative humidity, and CO 2 concentration are not sufficient to meet the breeding conditions, causing dairy barns to be exposed to low temperature and high humidity and poor air quality for a long time. Additionally, the cold stress evaluation findings suggest mild or moderate cold stress in dairy farming in the northern part of China. Although dairy cows are hardy ruminants, future adjustments to their winter housing and breeding environment are still needed (Galama et al., 2020).

Feasibility Analysis of Physiological Factors as an Index for Cold Stress Evaluation
The stress response is determined by a combination of animal and environmental factors, and the combination of the cow's physiological factors will provide a more realistic picture of its combined sensations (Fournel et al., 2017;Yan et al., 2019). Among the current cold stress evaluation methods, only WCT uses physiological indicators (skin temperature) as a modification of the cold stress index (Tew et al., 2002) but does not use physiological indicators as the influencing factors of cold stress in dairy cows. With the continuous development of sensor technology and animal behavior monitoring technology, various physiological characteristics of cows have been accurately obtained (Zhang et al., 2019). The RR is the most sensitive physiological factor to measure animal response to the external environment and can be considered an accurate factor to evaluate the degree of stress (Sun et al., 2021). Currently, researchers have proposed different methods to obtain RR. For example, Zhao et al. (2014) used the optical flow method to detect the RR of cattle because abdominal changes are obvious during lateral breathing, and the calculation accuracy was 95.68%. Yin et al. (2010) designed a cow behavior monitoring system with a wireless sensor node installed in the cow's neck to wirelessly acquire parameters such as the RR, body temperature, and locomotor behavior for the long-term observation of the cow's health status. The skin surface temperature of dairy cows is closely related to numerous physiological processes in the body and reflects health status (Khan et al., 2006). Peng et al. (2019) showed that the VST of cows is strongly correlated with rectal temperature. Additionally, the VST is more accessible than rectal temperature and can be used as a substitute for rectal temperature. Because of changes in the vasodilatory activity of the skin, the VST can be assessed by infrared thermography (Montanholi et al., 2008). This technology not only monitors skin temperature changes in cows (Suthar et al., 2012;Peng et al., 2019) but also measures udder health (Anagnostopoulos et al., 2021;Juan et al., 2021) and early lactation disease (Macmillan et al., 2019). Continuous measurement of physiological factors is not only an important tool for effective monitoring of animal health (Rutten et al., 2013;Hoffmann et al., 2019) but also provides data support for new methods for cold stress evaluation in dairy cows.
In the study, we analyzed the relationship between the RR and VST changes and milk yield in cows, and their highly significant correlations indicated that both could be relevant factors for cold stress evaluation. The average values of RR and VST for all cows were 23.39 ± 4.7 breaths/min and 23.90 ± 3.4°C, respectively. These data were slightly higher than the average RR (23 ± 5 breaths/min) and VST (21.21 ± 2.6°C) of cows in the winter derived from Chen (2014). The results of both studies indicated that under cold stress conditions in winter, cows must constantly reduce the number of breaths and body surface temperature to maintain thermal balance, which is a manifestation of the cow's adaptation to the environment.

Correlation Analysis of Behavioral Characteristics and Cold Stress
Because stress can adversely affect the lactation performance of dairy cows, the THI, BGHI, and ETI are used to classify the degree of cold/heat stress depending primarily on milk yield. Angrecka and Herbut (2015) found a strong correlation (r = 0.72-0.89 with P < 0.05) between milk yield and WCT in winter cows. To determine if there was a correlation between daily cold stress duration and milk yield, this study used a GLMM to statistically analyze the relationship between repeatedly measured daily cold stress duration and milk yield from the south and north sides. The results of this analysis showed that the correlation between the data was significant (all P < 0.001), and there were no significant differences between individual cows (all P > 0.05).
Behavioral changes in dairy cows are a response to external stress and are used to adjust imbalances in the body caused by physiological factors (Vaculiková et al., 2017;Shepley and Vasseur, 2021). When a cow is under prolonged cold stress, the cow will also act accordingly to maintain a constant body temperature (Redbo et al., 2001;Keren and Olson, 2006). Additionally, behavioral characteristics are also an outward indication of the health of the cow (Fournel et al., 2017). The behavioral responses of animals are used to regulate the housing environment and thus improve cow welfare (Redbo et al., 2001;Tucker et al., 2007;Angrecka and Herbut, 2015). The relationship between cold stress and behavioral characteristics is extensively investigated to facilitate timely detection of abnormal conditions in cows by breeders. Mader et al. (2010) analyzed the relationship between the CCI index and feed intake under cold conditions, and the correlation coefficient between the 2 was r ≥ 0.98, indicating a very strong correlation between the CCI index model and feed intake. In this study, significant correlations were also analyzed between daily cold stress duration and behavioral characteristics using GLMM to verify the validity of the cold stress evaluation method. The results of the statistical analysis showed that the feeding time and lying time increased and the active steps decreased when cows were under cold stress conditions, and these correlations were extremely significant (all P < 0.001). These findings were consistent with those of Redbo et al. (2001).
Feeding behavior is considered to be the main behavioral response of cows under cold stress (Méndez et al., 2020). When subjected to cold stress, cows automatically adjust by increasing DMI to increase their heat production and maintain body heat balance. Compared with the feeding time in the spring (5.79 ± 0.23 h/d) in the study, the feeding time on both the south and north sides increased by 11.4 and 25.9%, respectively. Additionally, Endres and Barberg (2007) and Ito et al. (2009) suggested that cow physical comfort is related to the lying and standing time, which represent the actual requirements of the cow to adapt to cold conditions. Gonyou et al. (1979) and Redbo et al. (1996) showed that cows spend more time lying down on cold days. The average lying time of individual cows in the present study was 11.9 ± 0.37 h/d (range 8.3 to 14.6 h/d), both higher than 11 h/d reported by von Keyserlingk et al. (2012) and 8 to 14 h/d reported by Haley et al. (2001). Conversely, cow active steps decreased with increasing cold stress duration. The reason is that standing increases the surface area of the cow's body exposed to the wind, increasing convective heat exchange with the environment (Vaculiková et al., 2017).
The cows selected for the present study performed well in terms of BCS (ranged from 2.75 to 3.5) and productive performance (Yukun et al., 2019). However, Marumo et al. (2022) showed that individual differences can lead to different responses of cows to external stresses. High producing cows are more susceptible to cold stress (Angrecka and Herbut, 2015). Tucker et al. (2007) also showed that lean cows changed their body temperature faster, took longer to feed and spent less time lying down than well-conditioned cows. Therefore, more extensive studies of the effects of cold stress on the physiological and behavioral characteristics of individual cows at different breeds and lactation stages are needed to ensure herd health in large-scale dairy barns and maximize productive performance.

CONCLUSIONS
In this study, a comprehensive evaluation index system was constructed using the AHP, and the GA was developed to optimize and establish index weights. Meanwhile, the evaluation index factors were further analyzed quantitatively and qualitatively, and a multilevel FCE method for cold stress in dairy cattle was proposed. The experimental results showed that cows on the south and north sides primarily experienced mild or moderate cold stress. Compared with the south side cows, which received more illumination, the north side cows experienced greater durations of mild and moderate cold stress (120.2 and 28.8 h, respectively). Sustained cold stress led to decreased milk production, increased feeding and lying time, and decreased standing time, especially for northern cows. The GLMM analysis showed that the daily milk yield and behavioral characteristics of dairy cows had a significant correlation with the daily cold stress duration (all P < 0.001). The results of this study reflect the cold stress experienced by cows in winter in northern China and provide a novel idea for evaluating the degree of cold stress.