Effect of BLG gene variants on milk-related traits in small ruminants: a meta-analysis

Abstract This study aimed to investigate the effect of A and B alleles of beta-lactoglobulin (BLG) gene on milk-related traits in small ruminants by using data reported in 26 potential articles. Records from 9439 genotyped sheep and goats were analysed by utilising four genetic models, including dominant (AA + AB vs BB), recessive (AA vs AB + BB), additive (AA vs BB), and co-dominant (AA + BB vs AB). The effect size of each genotype on daily milk yield, fat percentage, and protein percentage was measured as a standardised mean difference (SMD). The fixed-effects or random-effects model was used to analyse data based on the results of the heterogeneity test. The results showed that B allele of BLG gene positively, but slightly increased the fat percentage of milk in combined data of sheep and goat under dominant (SMD = 0.207, p = .018), recessive (SMD = 0.190, p = .000) and additive (SMD = 0.699, p = .003) genetic models. Also, a statistically significant effect of B allele on daily milk yield of goat was found under dominant (SMD = 0.208, p = .001) and co-dominant (SMD = −0.149, p = .002) models. There was no association between studied alleles and milk protein percentage in sheep and goat. Results of publication bias tests showed a lack of bias and sensitivity analyses suggested additive model as the best one and animals with BB genotype showed highest performance in milk fat production with a SMD equal to 0.699. Highlights Meta-analyses results indicated a small to medium association between BLG variants and daily milk yield in goats The B allele of BLG gene positively affect milk fat percentage in small ruminants The sensitivity analysis suggested additive genetic models for analysing data of fat percentage in small ruminants, based on the absence of outlier studies No association between BLG variants and milk protein was found


Introduction
Milk of small ruminants, such as sheep and goat, is an excellent source of nutritional content for people, especially in countries with limited water and pasture resources where cow rearing requires more cost and food resources compared to sheep and goat rearing. Hence, it has been the subject of considerable interest of many researchers to investigate the possible relationship between genetic polymorphisms of milk protein genes, milk yield and composition. It has been revealed that the milk of sheep and goat consists of two main milk proteins, namely caseins and whey proteins (Park et al. 2007). Whey proteins include a-lactalbumin (a-La) and b-lactoglobulin (b-Lg, BLG) and account for about 17-22% of the total proteins of milk (Mercier et al. 1978). BLG, the most abundant whey protein in the milk of ruminants, is secreted by mammary gland cells (Clark 1998). In sheep and goat, BLG gene is located on chromosome 3 and 11, respectively. Researchers have reported many allelic variants in promoter, coding and non-coding regions of this gene (Ballester et al. 2005;Sardina et al. 2012). Among them, a substitution of an adenine nucleotide by a guanine nucleotide at a specific position results in the replacement of Tyrosine amino acid with Histidine at the position 20 of the respected polypeptide (Kolde and Braunitzer 1983). Eventually, this substitution provides two different alleles named A and B that are known as main variants of BLG gene (Cardona et al. 2016;Jawasreh et al. 2019).
The impact of BLG variants on milk yield and composition of different breeds of sheep and goats has been widely investigated by researchers. It has been reported that these polymorphisms affect milk-related traits in sheep (C¸elik and € Ozdemir 2006;Dario et al. 2008;Padilla et al. 2018;Jawasreh et al. 2019) and goat (Kumar et al. 2006;Kahilo et al. 2014;El Hanafy et al. El Hanafy et al. El Hanafy et al. 2016;El-Shazly et al. 2017). However, some studies have suggested no effect of BLG variants on milk yield and composition traits in sheep (Giaccone et al. 2000;Kawecka and Radko 2011;Giambra et al. 2014;Triantaphyllopoulos et al. 2017) and goat (Kahilo et al. 2014;El-Shazly et al. 2017;Mahmood et al. 2016). Hence, this question could be raised whether the BLG variants influence milk-related traits or not.
It should be noted that some single studies have used small sample sizes to investigate the effects of BLG gene variants on milk yield and composition traits; therefore, making the conclusion from these studies may not be fully reliable. Moreover, the contradiction between the results of these studies is obvious. Hence, applying statistical techniques like meta-analysis could be a powerful tool to address the issues mentioned above by combining the data reported in all relevant studies and providing an extensive data set for analysing. Despite the increasing number of studies and data, the use of metaanalysis for investigating the effects of genetic variants on economic traits in animal science is limited (Mahmoudi et al. 2019;Mahmoudi et al. 2020). In this study, a meta-analysis on the effects of A and B alleles of BLG gene on the milk yield and composition of sheep and goat was performed by pooling data reported in eligible studies in the literature to evaluate the effects of BLG gene variants on milk related traits in small ruminants.

Literature search
To explore literature and selecting relevant studies to conduct a meta-analysis, we followed the Preferred Reporting Items for Systematic Reviews and Meta-Analyses (PRISMA) checklist criteria (Moher et al. 2009). All published journals and databases in different languages such as ScienceDirect, Wiley Online Library, PubMed, Springer, and Web of Science (in English), CNKI (in Chinese), Magiran (in Farsi), and eLibrary (in Russian) were explored independently by two authors up to Oct. 5 th 2020. The following search terms and their combinations were used to find appropriate studies: "b-Lg", "BLG", "LGB", "beta-lactoglobulin", "variant", "polymorphism", "sheep", "goat", "association", "relationship" and "milk". Furthermore, the reference list of found studies was scrutinised by authors to find out that no relevant article was missed in the meta-analysis.

Eligibility criteria and process of selection
We have considered articles as eligible for meta-analysis if they: (1) have investigated the effects of A and B alleles of BLG gene in the form of AA, AB and BB genotypes on the milk-related traits; (2) have reported sample size for AA, AB and BB genotypes of BLG gene; (3) have reported least-square means for each genotype; and (4) have reported standard deviation (SD), or standard error (SE) for the least-square mean of each genotype. The exclusion criteria were as follows: (1) duplicate articles; (2); articles in the form of an abstract; (3) articles with insufficient data; and (4) articles in the form of review.

Data extraction
To avoid conflicts between the first two researchers who selected studies for meta-analysis, a third researcher decided about inclusion or exclusion of identified studies. To extract data required for metaanalysis, two researchers independently checked all selected articles and added data in an Excel file for further calculations. Obtained data included the name of the first author of the article and the year of publication, species of animal (sheep or goat), the breed of sheep/goat, the sample size for AA, AB and BB genotypes, the least-square mean for each genotype and SD or SE reported for each least-square mean. In the case which SE was reported in the study, SD was calculated by using the following equation: where SE is the reported standard error of the mean for genotype, and N is the reported sample size for genotype. The equations suggested by the Cochrane Handbook for Systematic Reviews of Interventions (Higgins and Green 2011) were used to calculate pooled least-square means and SDs: where N 1 and N 2 are the sample sizes reported for two different genotypes, M 1 and M 2 are the leastsquare means and SD 1 and SD 2 are the standard deviations reported for the first and second genotypes, respectively. When there is an important study that is not reported the SD/SE, removing that study from the meta-analysis may affect the overall results of the meta-analysis. For this reason, we have used the imputation methods suggested by Furukawa et al. (2006) for imputing missing SD/SE in some studies.

Statistical analyses
In the present study, we have used the Review Manager v5.3 software package for analysing data. Collected data were divided into two groups based on species, including sheep (Merino, Lacaune, Awassi, Chios, East Friesian, Sarda, Heath, Valle del Belice, Bergschaf, Serra da Estrela, Leccese, Morkaraman, Churra, Massese, Zel, Polish Mountain, Altamurana and Manchega breeds) and goat (Sannen, Ardi, Hasbi, Harri, Barbari, Beetal, Jamunapari, Colombian tropical and African breeds). Furthermore, the combination of data from these two species was considered as a different group (total data). Because the p-value indicates the existence of an effect, but not the amount of that effect, providing effect sizes could be an appropriate way to clarify and interpret the results of the meta-analysis (McGough and Faraone 2009). In this meta-analysis, we have used the standardised mean difference (SMD) method to measure the effect sizes (Cohen 1988;Higgins and Green 2011). To interpret the effect sizes, Cohen defined the effect sizes of 0-0.2, 0.21-0.5, and more than 0.5 as small, medium, and large effects, respectively. We have employed four genetic models (Minelli et al. 2005;Thakkinstian et al. 2005;Lee 2015) including dominant (AA þ AB compared to BB), recessive (AA compared to AB þ BB), additive (AA compared to BB) and co-dominant (AA þ BB compared to AB) models for examining the relationship between BLG gene variants and milk-related traits in sheep and goat with 95% confidence interval (CI).
To quantify heterogeneity between studies, the I-squared (I 2 ) test with a range from 0 to 100% was conducted to ascertain that results of all studies are consistent (Lee 2015;Horita and Kaneko 2015). It has been suggested that I 2 less than 50% might not be important and could be considered as a low heterogeneity. In comparison, I 2 more than 50% indicates a high heterogeneity among studies and could have significant impacts on the results of meta-analysis (Vesterinen et al. 2014). The results of the heterogeneity test determine which model should be applied for analysing data. If there was heterogeneity, it indicated that the substantial diversity between studies and within-study variance should be considered. In this case, using a random-effects model could address the effects of diversity on the meta-analysis results. When there was no heterogeneity, it is assumed that only within-study variance exists and all studies are drawn from the same population; hence a fixed-effects model was used to calculate pooled effect sizes (DerSimonian and Laird 1986; Ades et al. 2005).
The Cook's distance measurement was used to estimate the influence of each data point and to find outlier studies (Viechtbauer and Cheung 2010). In the case of influential cases presence, we used the sensitivity analysis by removing one study at a time to prove that the overall results of the meta-analysis were not dependent on a single study. The bias of publication of articles included in the meta-analysis was evaluated by the Egger's test (Egger et al. 1997), and the funnel plots were checked to ensure that there was no publication bias.

Included studies and traits
The process of study selection is depicted as a PRISMA flow diagram in Figure 1. A total of 51 articles were found through journals and database searching along with their reference list of identified studies screening. From them, four articles were removed for being duplicated and three articles were in the form of an abstract; thus, they were also removed. In the next Figure 1. The inclusion and exclusion criteria as a PRISMA flow diagram. AA, AB and BB represent animals with two copy of A allele, one copy of A and one copy of B alleles, and two copy of B allele, respectively. step, 18 articles were excluded because: (1) they used different measurement such as breeding values instead of the least-square means; (2) they not focussed on BLG gene variants; and (3) they only reported the genotypic and allelic frequencies and not focussed on the relationship between BLG gene variants and milk-related traits. In total, 26 articles providing data of 9439 genotyped sheep and goats  (67) 6.24 ± 0.14 (100) 6.03 ± 0.25 (41) Yes Protein% 6.9 ± 0.3 (16) 6.78 ± 0.3 (16) 6.85 ± 0.2 (148) No a Least-square mean reported for each genotype. b Statistically significant results reported for the effect of BLG genotypes on the studied trait. c Daily milk yield (kg). d Not available. AA, AB and BB represent animals with two copy of A allele, one copy of A and one copy of B alleles, and two copy of B allele, respectively. Table 1) were considered to be eligible for the present meta-analysis. Some studies were investigated more than one breed; therefore, each breed was treated as a separate study in the present meta-analysis. The investigated traits included daily milk yield, milk fat percentage and milk protein percentage. The traits such as fat yield, protein yield, lactose percentage and solids-not-fat (SNF) percentage were ignored because the number of studies and sample sizes were not enough for conducting a meta-analysis. The characteristics of included studies investigated the effects of BLG gene variants on goat and sheep milk-related traits are presented in Tables 1 and 2, respectively.

Heterogeneity between studies
The I-squared (I 2 ) values obtained from the heterogeneity test under four genetic models and the number of studies used for performing meta-analyses of sheep, goat and their combined data are presented in Table  3. As can be seen, the I 2 for daily milk yield was less than 50%; hence, the fixed-effects model was used to analyse data of sheep, goat and their combination for this trait under all genetic models. Furthermore, studies that investigated the effects of BLG gene variants on fat and protein percentage in goats showed low heterogeneity, except for protein percentage under recessive model. In other cases, the random-effects model was used for analysing data.

Meta-analysis of BLG gene variants
The results of meta-analyses on the relationship between BLG gene variants and milk-related traits are presented in Table 4. Furthermore, the forest plots of different genetic models for studied traits are depicted in Figures 2-4.
The results showed that effect of genotype on daily milk yield was significant in goats under dominant (SMD ¼ 0.208, p ¼ .001, CI ¼ 0.080 to 0.336) and codominant (SMD ¼ À0.149, p ¼ .002, CI ¼ À0.241 to À0.256) models. However, the results indicated no association between genotype and daily milk yield in sheep and combined data, except for co-dominant .82 random a The fixed-effects and random-effects models were used, when I 2 was less than 50% and more than 50%, respectively. AA, AB and BB represent animals with two copy of A allele, one copy of A and one copy of B alleles, and two copy of B allele, respectively. AA, AB and BB represent animals with two copy of A allele, one copy of A and one copy of B alleles, and two copy of B allele, respectively. model (SMD ¼ À0.059, p ¼ .018, CI ¼ À0.108 to À0.010).
The association between genotypes of BLG gene and fat percentage was not observed in goats. However, results of meta-analyses in sheep showed a significant effect of genotype on fat percentage under dominant (SMD ¼ 0.104, p ¼ .003, CI ¼ 0.035 to 0.174), recessive (SMD ¼ 0.185, p ¼ .000, CI ¼ 0.115 to 0.254) and additive (SMD ¼ 0.227, p ¼ .000, CI ¼ 0.132 to 0.321) genetic models. In addition, meta-analyses of combined data of sheep and goat revealed impact of genotype on fat percentage under dominant (SMD ¼ 0.207, p ¼ .018, CI ¼ 0.035 to 0.378), recessive (SMD ¼ 0.190, p ¼ .000, CI ¼ 0.129 to 0.252) and additive (SMD ¼ 0.699, p ¼ .003, CI ¼ 0.230 to 1.169) genetic models. The relationship between BLG genotypes and protein percentage was found in neither goats nor sheep and their combination data under all genetic models.

Sensitivity analysis and publication bias
The results of data examination to detect influential studies in the meta-analyses are presented in Figure 5. It should be mentioned that only meta-analyses reporting significant results were tested for influential studies detection. The Cook's distance measurement introduced three studies including Corral et al. (2010), Giambra et al. (2014) and Jawasreh et al. (2019) as potent outliers. The sensitivity analysis was performed by omitting one of these studies at each pass to find out which one may have a very large impact on the overall outcome, and results indicated that meta-analysis outcome of the fat percentage in sheep and total data is sensitive to the study performed by Corral et al. (2010), and removing this study and re-analysing data resulted in non-significant outcome under all genetic models, except for additive model in total data. Besides, no publication bias was found for studies included in the present meta-analysis based on the results of Egger's test and shape of the funnel plots (data not shown).

Discussion
The results showed that BLG variants affected daily milk yield in goats under dominant and co-dominant models (p < .01), but there was no evidence of BLG variants influence under recessive and additive models (p > .05). From Figure 2, it can be concluded that B allele in the form of BB genotype slightly increases daily milk yield in goat. However, caution should be taken when interpreting these results, because the sample size available for goat is not large enough and this can affect the obtained results. There was no any indication of BLG variants' effects on daily milk yield in neither sheep nor the combined data of goat and sheep (p > .05). As can be seen from Table 4 and Figure 3, the fat percentage in the milk of sheep was strongly influenced by B allele of BLG gene under all genetic models (p < .01), except for co-dominant model (p > .05). The size of effect of B allele on fat percentage was low in dominant (SMD ¼ 0.104) and recessive (SMD ¼ 0.185) models, but high in additive model (SMD ¼ 0.227). It can be concluded that sheep carrying two B alleles produce milk with higher fat content in comparison with those that have only one B allele or two A alleles. Furthermore, the meta-analysis of combined data of sheep and goat confirmed the positive impact of B allele on fat percentage. In other words, despite the non-significant effect of B allele on fat percentage in the milk of goat, the overall results showed a relationship between B allele and fat percentage, especially when two alleles of B are carried together, which caused a high effect size (SMD ¼ 0.699). The non-significant finding in goats but the significant finding in combined sheep and goats is probably due to the overwhelming number of sheep in the metaanalysis. Hence, it could be a piece of strong evidence that the presence of B allele can increase the fat percentage in the milk of sheep. The meta-analysis results indicated that BLG gene variants did not affect (p > .05) protein percentage in the milk of sheep and goat, and meta-analyses of combined data confirmed these outcomes. As Figure 4 illustrates, the horizontal line touches the line of no effect, showing the nonsignificant effect of BLG gene alleles on milk protein percentage.
In a study performed on Awassi sheep, Jawasreh et al. (2019) reported that there was no difference between genotypes of BLG gene regarding daily milk yield. Furthermore, the BB genotypes of BLG gene statistically showed more protein percentage in comparison with AA and AB genotypes. They also reported that sheep having AA and BB genotypes produced more fat than those with AB genotype. In another study, Padilla et al. (2018) reported that Merino sheep carrying two copy of A allele produced more protein and fat percentage than those with one or none copy of A allele (AB and BB genotypes). Yousefi et al. (2013) reported that AB genotype of BLG gene showed more fat percentage than AA and BB genotypes in Zel sheep. However, they found no difference between abovementioned genotypes for protein percentage. In a study conducted on three Saudi goat breeds, El Hanafy et al. (2016) reported that goats having AA genotype of BLG gene produced more daily milk yield compared to those having AB and BB genotypes. Giambra et al. (2014) reported that identified variants on BLG gene have no significant effect on daily milk yield, protein percentage and fat percentage. From the results reported in previous studies, it is not clear which allele positively/negatively affects milk-related traits in sheep and goats; thus, the outcome of present meta-analysis can settle inconsistency between reported results.
As seen in Figure 5, all studies used in meta-analysis of BLG variants on daily milk yield under dominant and co-dominant models have the same influence and no influential study found. However, in the metaanalysis of fat percentage, three studies (showed in red circle) might be influential on the outcomes under dominant, recessive and additive genetic models in sheep as well as dominant and recessive models in sheep and goat combined data ( Figure 5(c-g)). However, Cook's distance cut-off (the black horizontal dotted line) did not pass by none of aforementioned studies under additive model in combined sheep and goat data ( Figure 5(h)), showing that the findings of meta-analysis using an additive model are not sensitive to these studies. Sensitivity analysis by omitting these studies once at a time indicated that only study performed by Corral et al. (2010) affected overall results. In other words, removing this study from meta-analysis resulted in non-significant effects of BLG variants on fat percentage under all genetic models in sheep and total data, except for additive model in total data.

Conclusions
In conclusion, the present meta-analysis showed that one copy of B allele of BLG gene may have a slight potential to increase fat percentage of milk in small ruminants; however, presence of tow copy of B allele results in a considerable increase in fat percentage of milk. Thus, the B allele of BLG gene could be used as a reliable marker for improving fat content in the milk of small ruminants in breeding programs, especially in sheep breeding. Nevertheless, the positive effects of BLG gene variants in enhancement of milk traits in goat needs more studies to be validated.

Ethical approval
Ethics approval was not required for this study.