**3. Processable Bulk Data Transfer (PBDT) tasks**

PBDT tasks require bulk transfer of processed data. Such data transfers are typical in multimedia systems and HEP experiments. For example in [1], 650MB of data was transferred on an average from a source to a set of sink nodes. High communication and computing times in PBDT tasks effectively amortizes the overhead of the LP-based algorithm used for optimization of the system performance. A PBDT task is characterized by the following three characteristics.

$$\mathbf{L} = \sum\_{\mathbf{i}=\mathbf{1}}^{\mathbf{k}} \mathbf{L}\_{\mathbf{i}} \tag{1}$$

Then for a PBDT task, the length of the required processed file is given by

$$\mathbf{L}^{\sim} = \boldsymbol{\Sigma}\_{\mathbf{l}=\mathbf{1}}^{\mathbf{k}} \boldsymbol{\varepsilon} \mathbf{L}\_{\mathbf{l}} \tag{2}$$

where ε*i* is a processing factor which is the ratio of the size of the processed partition and that of the original partition.

PBDT tasks are increasingly becoming important. They are used in various multimedia, high-energy physics and medical applications. The following section explains two of the practical examples of PBDT tasks.

### **3.1 Particle physics data grids**

Particle Physics Data Grids (PPDG) is a colloboratory project concerned with providing next-generation infrastructure for high-energy and nuclear physics experiments. One of the important requirements of PPDG is to deal with the enormous amount of data that is created during high-energy physics experiments that must be analyzed by large groups of specialists. Data storage, replication, job scheduling, resource management and security components of the Grid must be integrated for use by the physics collaborators. Processing these tasks require huge computing capabilities and fast communication capabilities. Grid computing is used for processing PPDG tasks that can be classified as a PBDT task.

### **3.2 Multimedia encoding**

52 Grid Computing – Technology and Applications, Widespread Coverage and New Horizons

GridBus workflow management effort [11] [23]. [23] has developed an architecture to specify and to schedule workflows under resource allocation constraints. Also, many of the data Grid projects that support distributed processing of remote data have proposed

PBDT tasks require bulk transfer of processed data. Such data transfers are typical in multimedia systems and HEP experiments. For example in [1], 650MB of data was transferred on an average from a source to a set of sink nodes. High communication and computing times in PBDT tasks effectively amortizes the overhead of the LP-based algorithm used for optimization of the system performance. A PBDT task is characterized by

1. The task involves large data transfer that has to be processed in some way before it can be used at the sink nodes. The large amount of data involved in the PBDT differentiates it from the compute intensive tasks where data usually consists of is only the parameters of the remote functions invoked. This implies that the data communication

3. The unprocessed raw file is such that it can be either processed as a whole or be divided into multiple partitions. If divided into partitions, each partition can be processed independently. The resultant processed partitions can later be combined to generate the required processed file. Consider a source file F, of size L. F can be partitioned into k

> L = ∑ L� �

 L~ = ∑ �L� �

where ε*i* is a processing factor which is the ratio of the size of the processed partition and

PBDT tasks are increasingly becoming important. They are used in various multimedia, high-energy physics and medical applications. The following section explains two of the

Particle Physics Data Grids (PPDG) is a colloboratory project concerned with providing next-generation infrastructure for high-energy and nuclear physics experiments. One of the important requirements of PPDG is to deal with the enormous amount of data that is created during high-energy physics experiments that must be analyzed by large groups of specialists. Data storage, replication, job scheduling, resource management and security components of the Grid must be integrated for use by the physics collaborators. Processing these tasks require huge computing capabilities and fast communication capabilities. Grid

computing is used for processing PPDG tasks that can be classified as a PBDT task.

��� (1)

��� (2)

workflow scheduling [11] [21].

the following three characteristics.

that of the original partition.

practical examples of PBDT tasks.

**3.1 Particle physics data grids**

**3. Processable Bulk Data Transfer (PBDT) tasks**

costs cannot be ignored while scheduling a PBDT task.

disjoint partitions, with data sizes of {L1, L2…. Lk}, such that

Then for a PBDT task, the length of the required processed file is given by

2. Cost of data processing is proportional to the length of the raw data file.

Multimedia encoding is required for applying a specific codec to a video [27]. Conventional methods use a single system for the conversion. The compression of the raw captured video data into an MPEG-1 or MPEG-2 data stream can take an enormous amount of time, which increases with higher quality conversions. Depending on the quality level of the video capture, the data required for a typical one hour tape can create over 10 GB of video data, which needs to be compressed to approximately 650 MB to fit on a VideoCD. The compression stage is CPU intensive, since it matches all parts of adjacent video frames looking for similar sub-pictures, and then creates an MPEG data stream encoding the frames. At higher quality levels, more data is initially captured and enhanced algorithms, which consume more time, are used. The compression process can take a day or more, depending on the quality level and the speed of the system being used. For commercial DVD quality, conversions are typically done by a service company that has developed higher quality conversion algorithms which may take considerable amount of time to execute. Grid technology is ideal for improving the process of video conversion.


The milk yield (kg/day) of every goat was determined during milking using a Lactocorder® device (Lactocorder, Balgach, Switzerland). This device collected a representative sample of 100 mL of milk at every milking of each animal for subsequent analysis. The macro-composition of milk (fat, protein, true protein, casein, whey protein, lactose, total solids, TS; non-fat total solids, NFTS; useful dry matter content, UDM, and ash; %) and urea content (mg/L) was determined by medium infrared spectroscopy equipment (MilkoScan™ FT2, Foss, Hillerød, Denmark). The SCC (103 cell/mL) was analysed by an electronic fluoro-optical method (DCC, DeLaval, Tumba, Sweden). Fat corrected milk yield (FCM) was calculated according to Gravert equation [33]: FCM (3.5%) = 0.433 × milk yield (kg/day) + 16.218 × fat milk yield (kg/day). Milk fatty acids were extracted by the Folch procedure, with some variations collected in Romeu-Nadal et al. [34] and were methylated following the Nudda et al. [35] method. The equipment, column, and FAME mix used for the identification of peaks of milk fatty acid profile were the same as for the diets. Atherogenicity index (AI) and thrombogenicity index (TI) were calculated according to Ulbricht and Southgate [36]. These indices provide important information because AI is related with the ability of lipids' adhesion to immunological and circulatory system cells and TI indicates the tendency to form clots in blood vessels [8]. Desaturase indices (DI) for C14:0, C16:0, and C18:0 were calculated according to Lock and Garnsworthy [37].

In order to assess the effect of the diets on goats' metabolism, blood samples were analysed. The same day as the milk sampling was performed, the fasting animals were bled and samples were collected for glucose, urea, and β-hydroxybutyrate (BHB) analysis. Blood samples were analysed with a glucose oxidase/peroxidase kit (Ref. 11503, Biosystems, Barcelona, Spain) for glucose (mg/dL), with a kinetic method (GN 10125, Gernon, Sant Joan Despí, Spain) for urea (mg/dL), and for the BHB (mmol/L), the Ranbut D-3-Hydroxybutyrate kit (RB 1007, Randox, Crumlin, UK) was used.

#### *2.4. Calculations and Statistical Analysis*

The SCC data were transformed into log10 scores before statistical analysis (LSCC).

BW, milk yield and macro-composition, SCC, and plasmatic profile data were performed using SAS GLIMMIX (SAS Institute Inc., Cary, NC, USA) with repeated measures, introducing the covariate of the data obtained in the pre-experimental sampling into the model and considering DIET, SAMPLING, and interaction DIET × SAMPLING as fixed effects, according to the following equation:

$$\mathcal{Y} = \mu + Di + S\mathfrak{j} + Dir\mathfrak{k}\mathfrak{j} + cov\mathcal{Y}O + Ak + \varepsilon\_r$$

where *Y* is the dependent variable, μ is the intercept, Di is the fixed effect of the diet (*i* = C, AP25, AP40, AP60, BB25, BB40, BB60), *Sj* is the fixed effect of sampling (*j* = 1, 2, 3), *DixSj* is the interaction of diet with sampling, *covY0* is the effect of the value of *Y* in sampling 0, *Ak* is the random effect of the animal, and *e* is the residual error. The covariance model with a lower value of the Akaike criterion (lower AIC and BIC) was used for each variable.

Milk mineral and fatty acid profile data were analysed using SAS GLM (SAS Institute Inc., Cary, NC, USA), introducing the covariate of the data obtained in the pre-experimental sampling into the model and considering DIET as a fixed effect. The level of acceptance for significance was 0.05.

#### **3. Results**

#### *3.1. Diet E*ff*ects on Body Weight and Feed Consumption*

Body weight is an indicator of the health status of the animal and optimising the inclusion of by-products involves assuring the proper health status of the goats. The treatments with the highest by-product inclusion showed a lower BW (40.2 and 38.7 kg in AP60 and BB60, respectively), while with the inclusion of 25% and 40%, no differences were observed compared to C (42.9 kg, Table 3). Sampling and interaction Treatment × Sampling also had a significant effect on BW as an increase (*p* < 0.001) was observed in sampling 2 in treatments with 40% of by-product (+1.9 and +2.4 kg in BB40 and AP40, respectively) and then in sampling 3, they descended again. Diets were offered in a similar amount but the goats in the different treatments showed different consumptions, with group C showing the highest (2.21 kg DM/day), whereas the lowest consumption was observed in groups BB40 (1.38 kg DM/day) and BB60 (1.27 kg DM/day) compared to the other treatments, which showed intermediate consumption (AP25: 1.52, AP40: 1.54, AP60: 1.57, and BB25: 1.65 kg DM/day).

 Body weight, milk yield, and composition and SCC, according to the effects considered.

C: Control diet; AP: Diet that includes artichoke plant silage; BB: Diet that includes broccoli by-product silage; 25, 40, and 60 inclusion level of by-product silage on dry matter basis %; SEM: Standard error mean; BW: Body weight; LSCC: Log10 somatic cell count; FCM: Fat corrected milk (3.5%); UDM: Useful dry matter content (fat + protein); TS: Total solids; NFTS: Non-fat total solids; abc: Least square means within a column having different superscripts differ significantly. \* *p* < 0.05; \*\* *p* < 0.01; \*\*\* *p* < 0.001.

#### *3.2. Milk Yield, Macro-Composition, and SCC*

A decrease in milk yield was observed as the percentage of inclusion of by-products increased (Table 3). C, AP25, and AP40 were the treatments with the highest milk daily yield (2.24, 2.15, and 2.14 kg/day, respectively; *p* < 0.001); BB60 was associated with the lowest yield (1.66 kg/day). A tendency to decrease FCM was also observed as the percentage of inclusion of the by-product in the diet increased. The highest yield was obtained in AP25, even without significant differences compared to C or other AP treatments; BB25 and BB40 did not show significant differences compared to C, AP40, and AP60, whereas BB60 showed the lowest value. The interaction among sampling and treatments was significant as the milk yield and FCM were only significantly reduced in AP25 and AP60 during the experiment, but remained stable in the rest of the treatments.

The diet had no significant effect on LSCC. An increase of + 0.28 Log cells/mL (*p* < 0.01) was observed in AP25 between samplings 2 and 3, so that sampling and interaction with treatment were significant.

As for the macro-composition of the milk shown in Table 3, the diet only had significant effects on fat (but also affected UDM and TS) and urea (Table 3). BB60 was the one with the highest fat value and T was the lowest. The significant interaction of the treatment with the sampling in fat, UDM, TS, whey protein, and lactose was due to specific increases or decreases in sampling 2 in AP40, which returned to similar values to the previous ones at sampling 3. Both the casein content of milk and NFTS were reduced in all treatments during the experiment (*p* < 0.001). The ash content increased 0.134 percentage units in AP25 at the end of the experiment, remaining stable in the rest of the treatments. Regarding the milk urea content, AP60 was the treatment that presented the highest level (641 mg/dL; *p* < 0.01) and BB60 the lowest (542 mg/dL).

#### *3.3. Milk Mineral Content*

Milk mineral profile is shown in Table 4. Only significant differences in the Mn concentration due to dietary treatment were observed, although of small magnitude. AP40 was the treatment that presented the highest level of Mn (0.233 mg/kg DM; *p* < 0.05), followed by BB25 (0.222 mg/kg DM), whereas BB40 was the treatment showing the lowest value (0.185 mg/kg DM). These differences between treatments are not considered biologically relevant because the greatest of them, which was between AP40 and BB40, was only 0.048 mg/kg DM.

**Table 4.** Milk mineral profile according to the effects considered.

C: Control diet; AP: Diet that includes artichoke plant silage; BB: Diet that includes broccoli by-product silage; 25, 40, and 60 inclusion level of by-product silage on dry matter basis %; SEM: Standard error mean. abc: Least square means within a column with different superscripts differ significantly. \* *p* < 0.05.

#### *3.4. Milk Fatty Acid Profile*

Regarding the milk fatty acid profile (Table 5), some significant variations were observed, although they were quantitatively limited. Regarding the content of vaccenic acid (C18:1t11), it was observed that this was higher (*p* < 0.001) in the AP treatments, without differences compared to C. There was a higher concentration of linoleic acid (C18:2n6) in AP60 (2.53%; *p* < 0.001); however, it was at C where a higher level of other C18: 2 isomers was observed. An increase (*p* < 0.001) of α-linolenic acid (C18:3n3) was observed as the level of AP inclusion in the ration was higher and AP60 presented a higher level (0.242%). AP treatments were also those with the highest content (*p* < 0.01) in rumenic acid (CLA c9, t11), although subsequently no significant differences were found in the sum of isomers of CLA (conjugated linoleic acid) between treatments, except of BB60, of which their content was the smallest of all. Table 6 shows that as the percentage of AP inclusion increased, so did the PUFA content, and AP60 was the treatment with the highest content (*p* < 0.001) compared to all the BB treatments, without differences from C or the rest of the AP treatments. AP60 presented the highest levels (*p* < 0.001) of n3 (0.275%) and n6 (2.79%) fatty acids, the latter without differences compared to C or the other AP treatments. It also achieved the lowest (*p* < 0.001) ratio n6/n3 obtained together with BB60 (10.3 and 12.3, respectively). Regarding the lipid quality indices related to human health (AI and TI), AP40 and AP60 were the ones with the lowest value (*p* < 0.001) and therefore, were healthier. Regarding the desaturation indices of the myristic (DI14), palmitic (DI16), and stearic (DI18) fatty acids, the differences found between treatments were of small magnitude. BB60 was the one with the highest value in DI14 and DI18 (0.012% and 2.08%, respectively; *p* < 0.001) and AP60 presented a higher value of DI16 (0.061%; *p* < 0.01).

**Table 5.** Fatty acid composition (g/100 g total fatty acids) measured in milk according to the effects considered.

C: Control diet; AP: Diet that includes artichoke plant silage; BB: Diet that includes broccoli by-product silage; 25, 40, and 60 inclusion level of by-product silage on dry matter basis %; SEM: Standard error mean; abc: Least square means within a column having different superscripts differ significantly. \* *p* < 0.05; \*\* *p* < 0.01; \*\*\* *p* < 0.001.

**Table 6.** Grouped fatty acids (g/100 g total fatty acids) and indices related to cardiovascular health and desaturation activity in milk according to the effects considered.

C: Control diet; AP: Diet that includes artichoke plant silage; BB: Diet that includes broccoli by-product silage; 25, 40, and 60 inclusion level of by-product silage on dry matter basis %; SEM: Standard error mean; SFA: Saturated fatty acids; MUFA: Monounsaturated fatty acids; PUFA: Polyunsaturated fatty acids; UFA: Unsaturated fatty acids (MUFA + PUFA); SCFA: Short chain fatty acids (C6:0 a C10:0); MCFA: Medium chain fatty acids (C11:0 a C17:0); LCFA: Long chain fatty acids (C18:0 a C24:0); AI: Atherogenic index; TI: Thrombogenic index; DI: Desaturation index; abc: Least square means within a column having different superscripts differ significantly. \* *p* < 0.05; \*\* *p* < 0.01; \*\*\* *p* < 0.001.

#### *3.5. Plasma Metabolic Profile*

Regarding the plasma metabolic profile (Table 7), it was observed that the greater the inclusion of BB in the diet, the higher the glucose level (49.5 and 50.0 mg/dL in BB40 and BB60; *p* < 0.001), although the differences were of small magnitude (42.5 mg/dL in BB25). Regarding urea, C and AP had a higher content (*p* < 0.001), while the BB treatments obtained lower levels and BB60 showed the lowest (33.2 mg/dL). The level of BHB was higher in treatments that included less by-product, such as AP25, AP40, and BB25, while it was lower in treatments that included more BB (0.299 and 0.304 mmol/L in BB40 and BB60, respectively; *p* < 0.001). There was significant interaction of treatment with sampling in the three variables due to the different behaviour throughout the experiment between treatments: Glucose increased (*p* < 0.001) with the progress of the experiment in all treatments except BB60; blood urea was reduced (*p* < 0.001) at sampling 2 in BB25 and BB40 and increased at sampling 3 in BB25, BB40, and BB60; BHB increased (*p* < 0.01) at the end of the experiment in BB25, BB60, and AP60, while in C, BB40, AP25, and AP40 remained stable.

**Table 7.** Plasmatic profile according to the effects considered.

C: Control diet; AP: Diet that includes artichoke plant silage; BB: Diet that includes broccoli by-product silage; 25, 40, and 60 inclusion level of by-product silage on dry matter basis %; SEM: Standard error mean; BHB: β-hydroxybutyrate; abc: Least square means within a column having different superscripts differ significantly. \*\* *p* < 0.01; \*\*\* *p* < 0.001.

#### **4. Discussion**

#### *4.1. Diet E*ff*ects on Body Weight and Feed Consumption*

One of the factors that affects the total volume of the diet and its consumption by livestock is the moisture content, as Jackson and Forbes [38] pointed out. This effect is especially important in the short term as herbivores are able to progressively modify the volume of the rumen to increase the speed of transit of the digesta [39], so in the long term, this effect would have less influence. In this experiment, carried out in the short term, diet C was the one presenting the highest DM content and feed consumption (2.21 kg DM/day). On the contrary, diets BB40 and BB60 contained a greater amount of water and were bulkier and presented less consumption. In addition, diets with silage showed higher concentrations of VFA and other substances resulting from fermentation. The presence of propionic acid in BB60 (4.79 g/kg DM), as well as a higher concentration of ammonia N in both BB40 and BB60, also occurred in treatments with lower consumption due to the depressing effect on feed consumption demonstrated by Baumont [40]. The feed consumption of the BB treatments was superior to those found by Meneses [41] (0.508 kg DM/day) in Murciano-Granadina castrated males, whose ration incorporated 55% of BB silage. All BW values were normal for the Murciano-Granadina breed [42,43]. The greatest reduction in BW was in BB60, as well as the greatest reduction in feed consumption (1.27 kg DM/day and 38.7 kg).

#### *4.2. Milk Yield, Macro-Composition, and SCC*

The treatments that presented a higher feed consumption were those that had a higher milk yield. The values obtained are similar to the yield obtained with the equation proposed by León et al. Ref. [44] for the modelling of the Murciano-Granadina lactation curve, which stands at 1.93 kg/day between the fourth and fifth lactation months, which is where the animals used in this experiment were located. The highest percentage of fat in BB60 (4.59%) was probably due to a concentration effect (being the treatment with the lowest yield) and its highest content in acetic acid (37.8 g/kg DM, triple the rest) in the diet, which is an extra-lipogenic nutrient precursor of fat synthesis. Van Knegsel et al. [45] observed similar effects in dairy cows when part of the corn in the diet was replaced by beet pulp. Due to a higher fat concentration in BB60, UDM and TS also reached the highest values in this treatment (8.03% and 12.9%, respectively), exceeding C by almost a percentage point. The urea level of all treatments was found to be within the optimal range for goats recommended by the Interprofessional Dairy Laboratory of Castilla-La Mancha (LILCAM), which is between 500 and 700 mg/L. The differences found in the milk urea content can be explained by the different levels of feed consumption of the treatments. BB60 presented less feed consumption, in particular refusing part of the offered BB, which probably induced lower total protein intake and lower levels of milk urea, as Jimeno et al. [46] noticed.

#### *4.3. Milk Mineral Content*

The macromineral values correspond to those found by Mellado and García [47] in goat crossings. The composition of the diet of animals largely determines the concentrations of minerals in milk [48]. As there were no large differences in the content of the different minerals in the diets, no significant differences were subsequently observed in the milk of the different treatments, which is important for the technological aptitude of the milk, given the relevance of Ca and P in the setting and development of the microstructure of cheese [49], the main destination of goat's milk. Only the Mn had a higher concentration in AP40 (0.233 mg/kg DM), although with such tight differences that they are not biologically relevant.

#### *4.4. Milk Fatty Acid Profile*

The milk of animals fed with AP60 had a higher content of n3 fatty acids, which caused a lower n6/n3 ratio, which is positive for the prevention of coronary and cardiovascular diseases [50]. On the other hand, C, AP25, AP40, and AP60, of which their diets had the highest levels of PUFA, were the treatments with milk richest in vaccenic, rumenic, and PUFA, as reported by Collomb et al. [51], who observed differences in the PUFA and vaccenic content in the milk of cows fed with high mountain pastures and in lowland plains because the plants that made up the mountain meadows had a higher concentration of PUFA.

Another factor that could influence the increase of PUFA in AP treatments was the slightly higher content of total polyphenols (TP) in the diet, although lower than that of BB60. However, the lower feed consumption of BB60 could mean that the total TP intake does not reach those of the AP treatments. Several studies have demonstrated the inhibitory action of dietary polyphenols on ruminal biohydrogenation of PUFA, without detrimental effects on milk yield and composition, due to interference with microbial flora [52–55]. These effects have also been observed in sheep with small amounts in the diet of by-products rich in TP [56,57]. Cappucci et al. Ref. [9] also observed how after increasing the TP content of the diet of Comisana sheep by including different levels of olive by-product, the concentration of linoleic (C18:2n6) and α-linolenic (C18:3n3) in milk was increased.

As a result of a lower milk content of C12:0, C14:0, C16:0, and C18:0, AP40 and AP60 had the lowest levels of AI and TI, so the milk of these animals would be of higher quality in terms of human health [42]. The values obtained from AI in all the treatments of this study are below those found by Molina-Alcaide et al. Ref. [42] in Murciano goats fed with conventional ration supplemented with feed blocks of olive by-products. The desaturation indices obtained in this experiment are similar to those provided by Baldin et al. Ref. [58] in a study in goats that received a dietary CLA supplement.

#### *4.5. Plasma Metabolic Profile*

Despite the differences found in the metabolic profile of the different treatments, glucose, urea, and BHB levels remained within the ranges considered optimal for goats [59], except for the urea value in BB60, which was slightly lower. As Friggens et al. [60] observed in goats' performance, the level of BHB was generally low and particularly in those treatments showing lower feed consumption (BB40 and BB60) because goats, as lactating animals, adapt their milk yield to the level of feed intake, as seen in Table 3. This reduces the metabolic load and allows them to maintain adequate body reserves turnover. Due to the strong relationship between plasma and milk urea content [61], the lower levels of blood urea were found in the same treatments with the lowest values of milk urea.

#### **5. Conclusions**

The findings of this study highlighted that a threshold level of AP or BB inclusion in dairy goat diets, without negative effects on milk yield, composition, mineral and fatty acid profile, as well as metabolic status of the animals, would be 40% of the dietary dry matter.

The inclusion of artichoke plant and broccoli by-product silages in high doses (60%) caused lower feed consumption and lower milk yield. Inclusion at 60% of AP and BB increased the milk TS, although not enough to compensate for the reduced yield, resulting in lower FCM in the case of BB60. No differences were found regarding the milk mineral profile. Inclusion of the artichoke plant silage in the animals' diet improved the milk lipid profile from the point of view of human health (AI, TI) compared to broccoli silage, due to a lower SFA content (C12:0, C14:0, and C16:0) and a higher concentration of PUFA, especially vaccenic acid (C18:1 trans11) and rumenic acid (CLA cis9, trans11), without any differences compared to the control treatment. Regarding sanitary status, the plasma metabolic profile in broccoli treatments reflects that goats ate grains and alfalfa, whereas broccoli was the last ingredient, impairing its consumption, especially at the high concentration (60%). In addition, the diets that included 60% of by-product silages caused a reduction in BW.

**Author Contributions:** Conceptualisation, G.R. and J.R.D.; Data curation, G.R., A.J.A.-B., and J.R.D.; Methodology, P.M., G.R., A.R., E.S., and J.R.D.; Resources, G.R.; Writing—original draft, P.M.; Writing—review & editing, G.R., A.S.A., C.A.S.-C., E.S., and J.R.D. All authors have read and agreed to the published version of the manuscript.

**Funding:** This research received funding from the Spanish Ministry of Economy, Industry and Competitiveness and the European Regional Development Fund, grant number AGL2015-64518-R (MINECO/FEDER, UE). Paula Monllor was funded by an FPU grant (Reference number: FPU14/06058) from the Spanish Ministry of Education.

**Acknowledgments:** Aprovertia S.L., a technology-based company of Miguel Hernández University of Elche, provided the facilities for manufacturing the silages.

**Conflicts of Interest:** The authors declare no conflict of interest. The funders had no role in the design of the study; in the collection, analyses, or interpretation of the data; in the writing of the manuscript, or in the decision to publish the results.

#### **References**

© 2020 by the authors. Licensee MDPI, Basel, Switzerland. This article is an open access article distributed under the terms and conditions of the Creative Commons Attribution (CC BY) license (http://creativecommons.org/licenses/by/4.0/).

MDPI St. Alban-Anlage 66 4052 Basel Switzerland Tel. +41 61 683 77 34 Fax +41 61 302 89 18 www.mdpi.com

*Foods* Editorial Office E-mail<EMAIL_ADDRESS>www.mdpi.com/journal/foods

MDPI St. Alban-Anlage 66 4052 Basel Switzerland

Tel: +41 61 683 77 34 Fax: +41 61 302 89 18

www.mdpi.com


*2.3. Extension, Relative Risk, Repair Cost, and Rainfall at Failure by Terrain Type*

First, by totaling the extensions of forest road segments by topographic type, we obtained the percentage of topographic types on an extension basis. The percentage of the landforms that were damaged was also calculated by adding up the extensions of the damaged segments by landform type. By comparing these data, the characteristics of failure in each landform type were examined in terms of the extension of failure.

Next, we calculated the relative risk of each landform type with respect to the other categories in order to examine the extent to which the susceptibility to failure differs by landform type. For example, the relative risk of the streamside category is calculated as the ratio of the probability of failure in the streamside category to the probability of failure in the other categories. Since the relative risk of each landform category represents the likelihood of failure relative to the other categories, it is also possible to compare the likelihood of failure among landforms. In order to clarify the relationship between topographic form categories and failure rates, a test of independence was conducted using the χ-square distribution. The null hypothesis of the χ-squared test is that the two variables, i.e., a given landform category and loss rate, are independent. If the null hypothesis is rejected, it suggests that there is some relationship between the two variables, i.e., the category may have an effect on the susceptibility to damage. Finally, we examined whether

the cost of repair and rainfall at the time of failure differed by terrain type. While repair costs and rainfall were recorded on a per-loss basis from the Forest Road Facility Failure Assessment Document, the topographic configuration differed on a per forest road segment basis. Therefore, depending on the length of the failure, a single failure may consist of segments of multiple terrain categories. In this case, the category with the highest number of segments was treated as the topographic category of the failure in question. The Kruskal– Wallis test was performed to determine whether the samples were drawn from the same population (or different populations with the same distribution) by comparing the median repair cost and rainfall at the time of failure for each topographic category group. If the null hypothesis is rejected, it suggests that not all samples may have been obtained from the same continuously distributed population. In this case, a further multiple comparison test is performed to determine whether there is a significant difference between any of the groups.

The rainfall at the time of failure recorded in the Forest Road Facility Failure Assessment Document was obtained from the nearest meteorological station (Automated Meteorological Data Acquisition System). These meteorological stations are located at intervals of about 17 km, and 45 of them are located in Nagano Prefecture (Figure 7) [31].

**Figure 7.** Distribution of meteorological stations established in Nagano Prefecture.

#### **3. Results**

#### *3.1. Percentage of Forest Roads Damaged by Segment Categories*

Figure 8 shows the percentage of all forest road segments and damaged forest road segments by terrain type category.

**Figure 8.** Percentage of all forest road segments and damaged forest road segments by terrain type category." Mixed" indicates segments that fall into both streamside and zero-order basin landforms.

The bar graph in the upper row of Figure 8 shows the percentage of forest road segment categories calculated for the 207 routes; 15% was streamside, 26% zero-order basin, 3% in both categories ("mixed"), 4% was stream crossing, and 53% was of other types.

The bar graph in the lower row of Figure 8 shows the percentage of forest road segment categories calculated for the failure segments; 42% was streamside, 18% zero-order basin, 8% in both categories, 7% was stream crossing, and 25% was of other types. By extension, 75% of the forest road failures that occurred in the target area occurred along streams or stream crossing or zero-order basins, indicating that forest road failures along streams are most common.

#### *3.2. Relative Risk of Forest Road Segment Categories*

The number of damaged and undamaged segments was recorded for each forest road segment category, and the relative risk to forest road segments in the "other" category was calculated (Figure 9). The χ2 test results showed that the ratio of the number of damaged to undamaged forest road segments in the "other" category was significantly higher than that in the "other" category for all categories. The relative risk for forest road segments in the "other" category was higher along streams and zero-order basins, in that order. The streamside was about 6.0 times more likely to be damaged, stream crossings were about 4.6 times more likely to be damaged, and zero-order basins were about 1.9 times more likely to be damaged than the other categories.

#### *3.3. Failure Characteristics by Forest Road Segment Category*

#### 3.3.1. Failure Repair Costs

Figure 10 shows the distribution of repair costs for the 287 damaged areas where repair costs were available (USD 1 = JPY 135.96 as of 7 March 2023). Only a small percentage of the failures were applicable to stream crossings (failures in which stream crossing segments were in the majority), so they were excluded from the analysis. The median and mean repair cost for all the damaged sites was approximately 2.33 million yen and 2.94 million yen, respectively (USD 16,946 and 21,584, respectively, as of 7 March 2023). The Kruskal–Wallis test was conducted for the three forest road segment categories, and the null hypothesis that all samples were derived from the same distribution with a 5% risk rate was not rejected. It is considered necessary to consider factors other than topographic morphology in analyzing the causes of damage that lead to high repair costs.

**Figure 9.** Cross tabulation table of forest road segment topography form and damage status. Rows correspond to the topographic form of the forest road segment, and columns correspond to the damage status.

**Figure 10.** Repair cost by risk category for each area damaged by failure. The line in each box represents the median of the sample; the upper and lower ends of each box represent the upper and lower quartiles, respectively.

#### 3.3.2. Amount of Rainfall

The distribution of the maximum 1 h rainfall and maximum 24 h rainfall for the rainfall events that triggered the forest road failures are shown in Figures 11 and 12, respectively. The maximum and minimum values for 1 h potential rainfall [32] and 24 h potential rainfall [33] at the target sites are plotted in the figures as reference values. The median maximum 1 h rainfall for all the damaged locations was 15.5 mm/h. The Kruskal– Wallis test was performed for the three forest road segment categories. However, the null hypothesis that all samples were derived from the same distribution with a 5% risk rate was not rejected. The quartile range of maximum 1 h rainfall at the damaged sites being located between the maximum and minimum values of 1-year potential rainfall suggests that forest road failures can occur even when the intensity of the maximum 1 h rainfall is not extreme. The median maximum 24 h rainfall intensity for all damaged sites was 175 mm/24 h. The Kruskal–Wallis test was conducted for the three forest road segment categories. The null hypothesis that all samples were derived from the same distribution with a 5% risk rate was rejected. Multiple comparisons showed that the mean maximum 24-h rainfall at the time of damage for the streamside category was significantly higher than the mean for the other categories at a 5% risk rate. The first quartile of the streamside category was higher than the other categories and exceeded the 30-year potential rainfall in some areas.

**Figure 11.** Distribution of maximum hourly rainfall in the damaged area. The line in each box represents the median of the sample; the upper and lower ends of each box represent the upper and lower quartiles, respectively.

**Figure 12.** Distribution of maximum 24 h rainfall in the damaged area. The line in each box represents the median of the sample; the upper and lower ends of each box represent the upper and lower quartiles, respectively.

#### **4. Discussion**

Streamside forest road segments accounted for only 15% of the total length of forest roads analyzed, but 42% of the road segments that were damaged. Furthermore, the relative risk of the streamside forest road segment was about 6.0 times higher than that of the other categories of forest road segments, indicating that it was the most likely terrain type to be damaged in this analysis. From the perspectives of both the length of failure and the relative risk, it is clear that the most important issue in the target area is the prevention of failure on the streamside forest road segments. The Forestry Agency supports the planning and improvement of trunk forest roads important for failure prevention and conducts projects to promote the strengthening of forest roads [34]. Countermeasures along stream segments will be particularly important for the resilience of Japan's forest roads against failures.

In previous studies evaluating the environmental impact of forest roads, distance to streams [35,36] and avoidance of stream crossings [37] have been used as one evaluation index. This study found that proximity to a stream is also an important evaluation indicator for forest road maintenance in terms of the susceptibility of forest roads to damage. Jing et al. [38] studied the spatial relationship between road and river networks in central China and noted that road and river networks were closely related spatially, with the density of high-standard roads increasing the closer they were to rivers. In Japan, there is a history of forest roads being located along streams because of the preference for timber haul-out routes using progressive gradients before the 1960s, when the driving performance of timber transport vehicles was low [39]. Especially when forest roads located along streams play a key role in the road network, they need to be subject to maintenance and improvement because of the significant disruption that would result from their destruction. The mean maximum 24-h rainfall at the time of the failure of a forest road segment along a stream was significantly higher than the mean for other geomorphic categories at a 5% risk rate. The frequency of extreme precipitation events per degree of warming also increased, and these trends were reported to be greater in Europe and Japan than in the United States and Australia [40]. Future changes in rainfall patterns may further increase the risk of damage, especially along streamside forest roads, and the relationship between failure and rainfall should continue to be investigated.

Even in terrain forms that were considered to be at high risk of failure in previous studies, the ease of failure in terms of relative risk differed by several fold. In studies of forest road evaluation based on multi-criteria evaluation, the AHP method has been used to determine the importance of different criteria [12,35–37]. Through interviews with experts and others, evaluations have been made, for example, that Criterion A is twice as important as Criterion B. However, these evaluations are based on the subjective opinions of experts. With regard to the weight of the criteria for failure, a more objective evaluation may be possible by accumulating statistical knowledge.

The topographic morphology of the forest road segments analyzed in this study could not explain the difference in failure repair costs. According to Watanabe et al. [41], who studied the past repair costs of 1504 forest roads in Nagano Prefecture, Japan, the standard deviation of the repair cost of one forest road failure amounts to about 26.2 million yen. Therefore, it is considered necessary to first investigate the characteristics of the failure in detail, especially for the failure with high repair costs.

In this study, we were able to obtain five years of inventory data on forest road failure by the prefectural government, which allowed us to statistically analyze the characteristics of forest road failure. In Japan, such inventory data are temporarily maintained by administrative departments but are discarded after a certain period of time because of a lack of familiarity with how to utilize them. As revealed in this study, some knowledge can be obtained only by collecting data on forest road failure on a scale of several hundred routes. Therefore, it is important to establish a system to compile a database of forest road failure inventories in the future. In this study, we focused on forest road failure that occurred more than 10 years ago, so we were unable to consider the presence, location, and condition of structures (drainage facilities, retaining walls, road surface pavement) that existed before the damage occurred. In Japan, it is rare to record information on such forest road structures in normal times because the value of their use is not well known. Since these structures are essential for a more accurate assessment of the risk of forest road failure, it is desirable to investigate their location and other information prior to the occurrence of a failure. Recently, studies have evaluated the degree of damage to failure-risk forest road surfaces by remote sensing [42–44]. Although it is labor-intensive to survey all forest roads immediately, it is also important to combine techniques to determine the functional status of forest road structures in a labor-saving manner, focusing on forest roads composed of segments at high risk of failure, as identified in this study and to build an inventory data of the structures.

There are various types of failure to forest roads, and it is necessary to construct an evaluation model of failure risk and prioritize countermeasures according to the cause and risk of failure. For example, in a forest road segment that runs alongside a stream, the nourishing force of the stream is considered to be the dominant failure factor, and reinforcement of the fill slope may be a more important failure risk reduction method than the placement of cross-drainage ditches. In this study, forest road segments were classified in terms of topography, and comparisons were made regarding the length of forest road failure, the relative probability of occurrence, repair cost, and induced rainfall intensity in each category. The results showed the importance of responding to forest road failures

along streams, which was revealed by the analysis based on a large number of data for a wide area of the prefecture.