1H HRMAS-NMR based metabolic fingerprints for discrimination of cheeses based on sensory qualities.

In this study, the 1H HRMAS-NMR (High-resolution Magic Angle Spinning-Nuclear Magnetic Resonance) spectra of 52 cheese samples obtained from the South Korean dairy farms were evaluated for their metabolic profiling and intensities associating with the sensory qualities. The NMR profiles displayed a broad range of compounds comprising amino acids, carbohydrates, organic acids, and phospholipids. Afterwards, the cheese samples were categorized into three groups (more likeness - G1, moderate likeness - G2, less likeness - G3), in relating to their sensory scores. The NMR data of the samples were later investigated through multivariate statistical tools to define the variations in metabolic fingerprints of every cheese sample and their intensities hailing in individual sensory groups. The unsupervised PCA employing all cheese samples unveiled the uniqueness in metabolite profiles of the brown and cheddar type cheeses (outliers). Moreover, Gouda and other types of cheeses displayed samples positioning in respective of their metabolite profiles. The pairwise comparison of sensory groups in the supervised models perceived better separation in OPLS-DA than PLS-DA. The corresponding VIP (PLS-DA) and loading (OPLS-DA) plots revealed amino acids and organic acids (lactate, citrate) as significant variables. The discrimination of G 1 Gouda type of cheeses against G 2 and G 3 was highly associated with their citrate levels. Further investigation using heatmaps displayed clear differentiation between each sensory group in terms of the levels of amino acids, lactate, citrate, phospholipids, and glycerol, conveying these variations are likely due to proteolytic and metabolic processes in cheese ripening. This study concluded that 1H HRMAS-NMR metabolite profile of the Korean cheeses is consistence with their sensory qualities. Further, the candidate metabolites identified in this study confers their potential application as biomarkers in cheese industries for faster and effective validation of sensory characteristics.

Introduction

Cheese is a highly popular fermented dairy product produced all over the world in numerous types for its specific flavour, aroma, and texture which are determined by the resources and production technologies precise to each region (Fox et al., 2017). Ripening is the most crucial step in cheese production, constituting a complex set of biochemical and microbiological reactions that give rise to the distinct sensory attributes in various cheese types. The ripening process and degree of ripening perform a crucial role in the quality improvement of cheese products with marginal cost and higher consumer acceptance (Santiago-López et al., 2018, Khattab et al., 2019). Additional factors, such as the raw materials, production process, and storage conditions, also significantly alter the cheese composition/quality (Ochi et al., 2012). Unfortunately, the domestic cheese industries still rely on traditional methods of grading and judging in quality control and product development. A precise study on the intact cheese composition is the utmost concern to obtain efficient quality control and describe the influence of individual variables. Some studies have focused on discriminating the cheese types based on specific compounds that affect distinct sensory attributes (Avsar et al., 2004, Van Leuven et al., 2008). These studies meant to be an alternative for the laborious, time-consuming, and inconsistent traditional methods.

The rising consumer awareness towards the quality and safety of food products obligates the industrial sectors to ascertain and affirm the identity, and description of the products (Drake and Delahunty, 2017). It is very complex in the cheese industry to manage cheese quality by grasping scientifically the potential changes in the production/ripening process and on the target compounds modifying the cheese quality (Ochi et al., 2013). Metabolite profiling holds promise as an extremely effective tool to explain the quality in terms of the ingredients, which would be more appropriate in the assessment of the cheese complex to ensure quality and reveal cheese production in general (Ochi et al., 2012, Scano et al., 2019a). Recently, the High-Resolution Magic Angle Spinning Nuclear Magnetic Resonance (HRMAS-NMR) spectroscopy aids as a powerful and versatile analytical tool for rapid and accurate identification of the overall metabolite profile in a single step (Mazzei and Piccolo, 2012). High reproducibility, non-selectiveness, relatively fast measurement, no sample extraction procedures, non-destructive and easier quantification of metabolites with lower cost under minimal analysis time are added advantages of HRMAS-NMR over other mass-spectrometry methods (Emwas et al., 2019). The metabolite profiles of cheese varieties of Parmigiano Reggiano (Shintu and Caldarelli, 2005), Emmental (Shintu and Caldarelli, 2006) and Mozzarella di Bufala Campana (Mazzei and Piccolo, 2012) are successfully analysed using HRMAS-NMR. In recent times, NMR profiles of cheeses combined with multivariate statistical approaches implemented successfully to discriminate their geographical origin (Consonni and Cagliani, 2008), age of ripening periods (Shintu and Caldarelli, 2005), brining time (Ruyssen et al., 2013) and adjunct cultures (Rodrigues et al., 2011, Piras et al., 2013).

The continuous demand for domestic cheese consumption in both the consumer market and food processing industry increased the proportion of local dairy farms for cheese production in South Korea (Kandasamy et al., 2019). To date, no literature published regarding the metabolomics studies on Korean cheese types and in particular, no NMR profile data available or remains scarce. In present study, the metabolite profiling of different cheese types of South Korea was determined through 1H HRMAS-NMR analysis and evaluated for their discrimination based on their sensory qualities, with the ultimate aim of acquiring a powerful tool for rapidly tracing the quality of cheeses to consumer perception. This study focused on (i) to reveal the differences of the metabolite profiles in various cheese types through 1H HRMAS-NMR analysis; (ii) to evaluate the variations in sensory qualities of cheeses by coupling the NMR data with multivariate statistical analysis; and (iii) to classify the potent metabolites that alter the acceptable characteristics in cheeses.

Materials and methods

Cheese samples

Fifty-two cheese samples, including 11 types (Gouda [34], Camembert [5], Berg [3], Appenzeller [2], Cheddar [2], Brie [1], Brown [1], Emmental [1], Frill [1], Quark [1], and Reblochon [1]), collected from the local dairy farms of seven different provinces in South Korea were the basis of this study. The numbers in parenthesis represent the number of samples in each cheese type, and their ages ranged from 3 to 6 months. The complete details of the number of samples in each cheese type and their geographical origin are available in Table 1.

Table 1

Details about the cheese types and location of all the samples used in study.

Sample Id.	Cheese type	Location	Province	Total
1	Gouda	Paju-si	Gyeonggi-do	18
2	Berg	Pocheon-si
4	Gouda	Goyang-si
6	Gouda	Goyang-si
17	Gouda	Yeoncheon-gun
20	Gouda	Goyang-si
23	Berg	Paju-si
26	Gouda	Pocheon-si
29	Gouda	Paju-si
30	Gouda	Goyang-si
34	Appenzeller	Paju-si
35	Gouda	Yeoncheon-gun
38	Gouda	Goyang-si
39	Gouda	Paju-si
42	Gouda	Paju-si
43	Gouda - black pepper	Paju-si
50	Gouda	Goyang-si
51	Gouda	Yeoncheon-gun
3	Gouda	Cheonan-si	Chungcheongnam-do	9
22	Gouda	Cheonan-si
28	Camembert	Cheonan-si
36	Cheddar	Taean-gun
37	Gouda	Cheonan-si
40	Camembert	Cheonan-si
46	Gouda	Cheonan-si
48	Camembert	Cheonan-si
12	Gouda	Cheonan-si
9	Gouda	Yeonggwang-gun	Jeollanam-do	13
11	Gouda	Jangheung-gun,
52	Brown	Yeonggwang-gun
13	Gouda	Yeongam-gun
14	Camembert	Yeongam-gun
18	Berg	Yeongam-gun
19	Gouda	Yeonggwang-gun
24	Emmental	Yeonggwang-gun
25	Pepper camembert	Yeongam-gun
32	Gouda	Yeongam-gun
33	Quark	Yeonggwang-gun
41	Gouda	Hampyeong-gun
44	Gouda	Yeongam-gun
49	Brie	Yeongam-gun
5	Gouda	Hamyang-gun	Gyeongsangnam-do	5
16	Frill	Hamyang-gun
21	Reblochon	Hamyang-gun
31	Bongson Cress (Gouda)	Hanan-gun
47	Gouda	Hanan-gun
10	Cheddar	Gimje-si	Jeollabuk-do	4
15	Gouda	Jeongeup-si
27	Appenzeller	Gimje-si
45	Gouda	Gimje-si
7	Gouda	Cheorwon-gun	Gangwon-do	1
8	Gouda	Cheongju-si	Chungcheongbuk-do	1

Sensory scoring of cheeses

Sensory evaluation was conducted in cheeses using six expert personnel from different universities, research institutes and hotels in South Korea. The experts are well qualified in carrying sensory analysis at a high level of discrimination, sensitivity and consistency in measurement. The cheese samples were temporarily stored in a refrigerator and placed at room temperature one hour before tasting. The samples were evaluated by the experts for sensory scores (Table 2) according to flavours (acidic, bitter, barny, cooked, feed, fruity, fermented, lipase/rancid, mouldy, weedy, yeasty), body and textures (curdy, gassy, mechanical holes, pasty, sponge, weak, abnormal/ heterofermentative), and appearance and colours (cracks in the rind, surface mould, rough surface, soiled surface, soft spots, huffed). De-ionised water and unsalted crackers were provided for palate cleaning between each sample. The score data of all cheese samples were managed and analysed using Microsoft Excel spreadsheet.

Table 2

Sensory score chart for evaluation of cheese samples.

Sensory properties	Score
Flavour	50
Body & Texture	30
Appearance & Colour	20
Total	100

Sample preparation for NMR analysis

3-(trimethylsilyl) propionic-2,2,3,3-d4 acid sodium salt (TSP-d4) and deuterium oxide (D2O) were obtained from Sigma-Aldrich Co (St. Louis, MO, USA). For NMR analysis, each cheese sample weighing 20 mg added into a 4 mm NMR nanotube followed by addition of 20 μL of 2 mM TSP-d4 solution in D2O.

NMR data acquisition and processing

The NMR spectral measurements were obtained on a 600 MHz Agilent HR-MAS (High Resolution-Magic Angle Spinning) NMR spectrometer (Agilent Technologies, Palo Alto, CA, USA) running at a 1H frequency of 599.93 MHz, outfitted with a 4-mm gHX NanoProbe by a spinning rate of 2000 Hz. The presat Carr–Purcell–Meiboom–Gill (CPMG) pulse sequence employed to control the residual water signal and high-molecular-mass compounds. Throughout the experiments, the sample temperature constantly maintained at 25 °C (298.15 K). For each sample, 128 transients acquired with a relaxation delay of 3.0 s, an acquisition time of 3.0 s, the proton 90° pulse (7.075 μs) which resulted in a total acquisition period of 13 min 9 sec for each sample. The entire spectral data were Fourier-transformed after multiplying the free induction decays by an exponential weighting function corresponding to a line broadening of 0.5 Hz.

Targeted metabolite profiling

Identification of metabolites was performed using the Chenomx Profiler, a module of Chenomx NMR Suite 7.1 professional edition (Chenomx Inc., Edmonton, Canada). While, the quantification of metabolites were done using the Chenomx 600 MHz library database from Chenomx NMR Suite 7.1, in which the concentration of individual metabolites was determined through the signal corresponding to a known concentration of reference (in this case TSP-d4). Each spectrum was converted to the frequency domain, phased and baseline corrected, and then the TSP-d4 singlet peak was calibrated to 0.00 ppm.

Multivariate statistical analysis

To analyse the variations in the metabolite profile of cheese samples and within the sensory groups, multivariate statistical analysis methods viz., PCA (Principal Component Analysis), PLS-DA (Partial Least-Squares Discriminant Analysis) and OPLS-DA (Orthogonal PLS-DA) analysis were conducted in the webserver MetaboAnalyst 4.0 (Chong et al., 2019). Initially, the 1H NMR data were log-transformed, normalized and Pareto scaling (divide the mean-centred data by the square root of the standard deviation for the variable) applied to ensure an unbiased variability between the samples. Pareto scaling is generally applied when a large dynamic variation exists in the data set. Compared to unscaling, Pareto scaling gives greater weight to the variables in the generated models by increasing the contribution of metabolites with a lower concentration (Craig et al., 2006).

Among the multivariate models, PCA, an unsupervised dimensionality reduction method was executed initially to detect the possible outliers and the intrinsic clustering among the groups, without any bias. Subsequently, the supervised models PLS-DA and OPLS-DA were performed to obtain a maximum separation between classes. PLS-DA is a linear classification model predicted based on accuracy, multiple correlation coefficient (R2), and cross-validated R2 (Q2) in 10-fold cross-validation. R2 represents a quantitative measure of the total variation in the data indicating the goodness of fit and Q2 the goodness of predictability (Worley and Powers, 2013). The major distinctive features were ranked based on the scores (≥1) of the variable importance in projection (VIP) of each variable in the PLS-DA model. While OPLS-DA is a modified PLS-DA method, that explains the variation between and within the groups, using distinct predictive and orthogonal components. The S-plot (scatter plot) generated in this model, visualize the covariance (p) and the correlation structure p(corr) among the X-variables and the predictive score t[1]. This model was estimated using R2X (explained variation in X), R2Y (explained variation in Y), and Q2Y (predicted variation in Y) using cross-validation. The values of these parameters are between 0 and 1; the nearer they approach 1, they can be better predicted or explained.

Further, using the normalized data, an HCA heat map was produced to describe the classification ability and the concentration of each metabolite in all the cheeses. In addition, heat maps based on the sensory groups were constructed to determine the metabolite profile difference within the groups.

Results and discussion

Sensory evaluation of cheeses

Sensory evaluation of cheese is of utmost important for ascertaining the virtual merits in the cheese-making practices and the impact of constituents on the precise sensory traits of cheese. Besides, sensory evaluation is highly essential to govern the effect of sensory characteristics of the intake attributes of cheese and consumer satisfaction (Drake and Delahunty, 2017). In this study, the sensory evaluation was conducted through a group of cheese experts to calculate the score of acceptability for attributes viz., flavour, body and texture, and appearance and colour, separately for every cheese sample as described in Section 2.2. The reliable quantitative measurement on the evaluation of sensory characteristics is vital in the food industry to drive the fluctuating needs towards both consumers and the market (Khattab et al., 2019). Based on an overall score assessment, the cheese samples were classified under three different sensory groups (G1, G2, G3) as outlined in Table 3. The cheese samples scoring 507–457 (G1) was rated as more likeness, while 456–406 (G2) as moderate likeness and <406 (G3) as less likeness. Among the 11 cheese types evaluated, the samples of Gouda were predominant and distributed in all the groups. Interestingly, all Gouda cheese samples of Chungcheongnam-do province is rated as G1. On the contrary, other cheese samples collected in a single province found disseminated in all three groups. These results suggest that the sensory attributes of Korean cheese samples are not greatly affected by their geographical resources.

Table 3

Details about sensory scoring and group classification of cheese samples used in multivariate analysis.

Group & Sensory Score	Cheese variety	Location	Samples (nos)
I(507–457)	Appenzeller	Gyeonggi-do	1	13
	Frill	Gyeongsangnam-do	1
	Gouda	Chungcheongnam-do	5
		Jeollanam-do	2
		Jeollabuk-do	2
		Gyeongsangnam-do	1
	Quark	Jeollanam-do	1
II(456–406)	Appenzeller	Jeollabuk-do	1	23
	Brie	Jeollanam-do	1
	Berg	Gyeonggi-do	2
	Camembert	Chungcheongnam-do	1
		Jeollanam-do	1
	Camembert (Black Pepper)	Jeollanam-do	1
	Cheddar	Chungcheongnam-do	1
	Gouda	Gyeonggi-do	8
		Jeollanam-do	3
		Gangwon-do	1
	Gouda (Black Pepper)	Gyeonggi-do	1
	Gouda (Bongson Cress)	Gyeongsangnam-do	1
	Reblochon	Gyeongsangnam-do	1
III<406	Berg	Jeollanam-do	1	16
	Brown	Jeollanam-do	1
	Camembert	Chungcheongnam-do	2
	Cheddar	Jeollabuk-do	1
	Emmental	Jeollanam-do	1
	Gouda	Gyeonggi-do	6
		Jeollanam-do	2
		Gyeongsangnam-do	1
		Chungcheongbuk-do	1

1H HRMAS-NMR analysis and metabolic assignments of cheese samples

In cheese ripening process, a cascade of microbiological and biochemical events (glycolysis, lipolysis and proteolysis) preceded by the starter and adjunct cultures gives rise to several bioactive compounds that control the organoleptic characteristics (Santiago-López et al., 2018, Scano et al., 2019a). Metabolite profiling of the intact cheese using HRMAS-NMR allows the direct documentation of metabolites without prior extraction processes (Mazzei and Piccolo, 2012). A representative 1H HRMAS-NMR spectrum obtained from the cheese samples employed in this study (Fig. 1), with a total of 29 identified metabolites. The overall information concerning the proton NMR assignments, chemical shifts and peak multiplicity is available in Table 4. The NMR spectra identified 16 amino acids (alanine, asparagine, aspartate, glutamate, glutamine, glycine, isoleucine, leucine, lysine, methionine, phenylalanine, proline, serine, threonine, tyrosine and valine); 5 organic acids (acetate, citrate, lactate, pyruvate and succinate; 2 sugars (galactose and lactose); 3 amino acid derivatives (choline, creatine and tyramine); 2 phospholipids (O-phosphocholine and sn-Glycero-3-phosphocholine); and an alcohol of glycerol. Similar, metabolite profiles are previously compiled in Parmigiano Reggiano (Consonni and Cagliani et al., 2008) and Fossa pit (Scano et al., 2019b) cheeses. Of all the metabolites, proportion of lactate and aminoacids are predominant in all the cheese samples. Ochi et al., (2013) stated that the levels of aminoacids and lactate aids as strong marker candidates in cheese ripening process. Similar 1H NMR cheeses profile with lactate as the most abundant metabolite is stated earlier by Ruyssen et al. (2013).

Fig. 1

Representative 600 MHz 1H HRMAS NMR spectrum of complete (a) and expanded regions (b) of the metabolites in cheeses. The assignments to the functional molecules are reported in Table 4.

Table 4

1H HRMAS-NMR assignments, identified metabolites and chemical shifts (multiplicitya) of distinguishable peaks obtained in D2O of cheese samples.

Assigned No.	Metabolite	Chemical shift (ppm)
1.	Acetate	1.9 (s)
2.	Alanine	1.5(d), 3.8(m)
3.	Asparagine	2.8(q), 2.9(q), 4.0(q), 6.9(s)
4.	Aspartate	2.7(q), 2.8(q), 3.9(q)
5.	Choline	3.2(s), 3.5(t), 4.1(m)
6.	Citrate	2.5(d), 2.7(d)
7.	Creatine	3.0(s), 3.9(s)
8.	Galactose	3.5(m), 3.6(m), 3.7(m), 3.8(m), 3.9(m), 4.0(m), 4.1(m), 4.6(d), 5.3(d)
9.	Glutamate	2.1(m), 2.4(m), 3.8(m)
10.	Glutamine	2.1(m), 2.4(m), 3.8(t), 6.9(s)
11.	Glycerol	3.5(q), 3.69(m), 3.8(m)
12.	Glycine	3.5(s)
13.	Isoleucine	0.9(t), 1.0(d), 1.2(m), 1.5(m), 2.0(m), 3.7(d)
14.	Lactate	1.3(d), 4.1(q)
15.	Lactose	3.3(t), 3.5–4.4(m), 4.7(d), 5.2(d)
16.	Leucine	0.9(t), 1.7(m), 3.7(m)
17.	Lysine	1.4(m), 1.5(m), 1.7(m), 1.9(m), 3.0(t), 3.9(t)
18.	Methionine	2.1(s), 2.2(m), 2.6(t), 3.8(m)
19.	O-Phosphocholine	3.2(s), 3.6(m), 4.2(m)
20.	Phenylalanine	3.1(q), 3.3(q), 4.0(m), 7.3–7.4(m)
21.	Proline	2.0(m), 2.3(m), 3.3 (m), 3.4(m), 4.1(m)
22.	Pyruvate	2.4(s)
23.	Serine	3.8(m), 3.9(m), 4.0(m)
24.	Succinate	2.5(s)
25.	Threonine	1.3(d), 3.6(d), 4.3(m)
26.	Tyramine	2.9(t), 3.2(t), 6.9(d), 7.2(d)
27.	Tyrosine	3.0(m), 3.2(m), 3.9(m), 6.9(d), 7.2(d)
28.	Valine	1.0(q), 2.3(m), 3.6(d)
29.	sn-Glycero-3-phosphocholine	3.2(s), 3.6(m), 3.9(m), 4.0(m), 4.3(m)

Letters in parenthesis denote the peak multiplicities: s - singlet; m - multiplet; d - doublet; t – triplet and q – quartet.

Multivariate analysis

The relevance of HRMAS-NMR metabolite profiles in the differentiation of ripening periods, geographical origin and quality assessment has been demonstrated earlier in Parmigiano Reggiano (Shintu and Caldarelli, 2005), Emmental (Shintu and Caldarelli, 2006) and Mozzarella di Bufala Campana (Mazzei and Piccolo, 2012) cheeses. In the current study, we combined 1H HR-MAS NMR with multivariate statistical methods to examine whether the metabolite fingerprints of cheese samples are suitably unique to detect metabolic markers in the sensory group discrimination. To our knowledge, this is the first study to deal with the metabolite profiles and intensities of the Korean cheeses using 1H HRMAS-NMR and to discriminate them based on their sensory qualities.

PCA analysis for cheese samples

The 1H HRMAS-NMR metabolite data matrix obtained for the cheese samples was assessed initially by Principal component analysis (PCA) to reduce the dimensionality of the data, visualize samples grouping and outliers. This untargeted PCA aids in the detection of any inherent sample clustering devoid of any bias since it does not require any information on the data sets (Worley and Powers, 2013). The initial PCA assessed the metabolite profiles of all the 52 cheese samples, divided into two classes. The class 1 and 2 represents Gouda type and other types of cheeses, respectively. The PCA score plot (Fig. 2A) shows large overlapping and tight clustering among cheese samples in the centre, with total variance accounting 45.9% as the first principal component (PC1) and 14.7% as the second principal component (PC2). Nevertheless, one strong (52 - brown cheese) and one weak (36 - cheddar cheese) outliers remain outside and on the Hotelling T2 ellipse (95% confidence region) respectively, are identified for uniqueness in their metabolite profiles. The brown cheese shows positive contributions on both the axis, while cheddar cheese shows positive in PC1 and negative in the PC2 axis. The loading plots (Fig. 2B) shows creatine, galactose, glycerol, lactose, O-phosphocholine and threonine were positive towards both the PCs explaining their higher levels in brown cheese. Likewise, choline, pyruvate and citrate are positive in PC1 and negative in PC2 indicating these variables define the composition of cheddar cheese.

Fig. 2

Score (left) and loading (right) plots of Principal Component Analysis (PCA) showing the metabolic pattern for all the cheese samples (A, B; 1-Gouda type, 2-other types of cheeses), only Gouda type (C, E) and remaining types of cheese samples (D, F). In score plots of C and E, 1, 2, 3 represent the groups classified based on sensory evaluation.

Subsequent PCA analysis after excluding the outliers was challenging to elucidate since the samples in the class of other types of cheeses were detected as outliers consistently, may be due to the dominance of Gouda type cheeses. Hence, we ran a PCA analysis separately for the classes Gouda type and other types of cheeses, to study their relevant discrimination based on sensory qualities (see Section 3.1).

In the PCA score plots of both Gouda type (Fig. 2C) and other types (Fig. 2D) of cheeses, all samples were lying inside the Hotelling T2 ellipse. These plots demonstrated a clear view that each cheese sample is plotted based on the NMR metabolite composition, rather than on geographical location. The score plot for the Gouda type cheeses exhibits a total variance accounting of 39.3% in PC1 and 16.1% in PC2 (Fig. 2C). The sensory G1 Gouda samples are noticed to be closely clustering on the negative side of PC1 and positive on PC2. The samples closely clustering together signify similar metabolite composition, while more distance indicates greater variation. No clear separation detected among the sensory G2 and G3 samples. In the loading plot (Fig. 2E), the amino acids mainly negative correlated in PC1 correlates the higher proportion of amino acids in G 1 cheese samples. The distantly found variable citrate contributes negative in PC2 and positive in PC1, outlining its significance in the separation between sensory groups of Gouda type cheeses. In Gouda cheese production, citrate acts as a main source for diacetyl (2,3-butanedione) production by the citrate metabolising bacteria employed as starter cultures (Van den Berg et al., 2004). The phospholipids and lactate variables highly influence the composition of the Gouda-Bongsen Cress (31) that lies near the boundary of the Hotelling T2 ellipse.

Likewise, in PCA for class other types of cheeses (except Gouda type), the score plot exhibits a total variance accounting of 34.1% for PC1 and 16.2% for PC2 (Fig. 2D). In the score plot, the Camembert cheese (40) remains distant from other samples with a positive contribution to both the PCs. The variable tyramine located in far and having a positive contribution on both the PCs in the loading plot (Fig. 2F) is bound for the uniqueness in Camembert cheese.

PLS-DA of cheese samples

Next, we employed a supervised PLS-DA approach, which maximizes the separation between the groups by presenting the highest predicted indicator variable. This model on Gouda type cheeses improved better separation between the groups with a total variance of 39% and 7.6% in PLS components 1 and 2, respectively. Compared with PCA, a clear separation of sensory G1, and a certain degree of separation between G2 and III samples are displayed in this model (Fig. 3A). The location of G1 cheese samples along the direction of PLS component 2 is attributed with their higher content of organic acids (acetate, succinate, aspartate) and amino acids together with lower levels of glycerol, O-phosphocholine, choline, citrate, lactate and sn-Glycero-3-phosphocholine (Fig. 3B). The statistical robustness of the model evaluated by 10-fold cross-validation shows the determination coefficient (R2) and cross-validation determination coefficient (Q2) was 0.665 and 0.494, respectively. Further cross-validation of the model done through 20 permutation tests gives a significant p-value of <0.05. The Variable importance projection (VIP) plot derived in this model provides a holistic view of significant variables ranked from high to low VIP scores. The score of each variable is calculated as the weighted sum of the squared correlation between the initial variables and PLS-DA components. Variables with VIP scores around 1.0 works well in selection of significant variables, although scores >1.0 are most relevant in explaining the model (Worley and Powers, 2013). Fig. 3C shows the metabolites such as citrate, isoleucine, lysine, threonine, valine, succinate, proline, O-phosphocholine, glycine, glutamate, lactate, choline, asparagine and methionine as the main discriminants in the PLS-DA model. Among them, citrate yields the largest VIP score attributing to its highest contribution in sensory group discrimination.

Fig. 3

PLS-DA score, loading and VIP plots derived from the 1H NMR for all the cheese samples in Gouda type (A-C) and samples except Gouda type (D-F). In the score plots (A, D), 1, 2, 3 represent the groups classified based on sensory evaluation.

Similar PLS-DA approach for samples of other types of cheese accounted for a total variance of 25.8% in PLS 1 and 17.0% in PLS 2 components (Fig. 3D) with appreciable group discrimination between G1 and G3 samples. Higher citrate, lactate and tyramine levels together with lower amino acids in G3 samples and the vice-versa of this composition in G1 define the separation between the groups (Fig. 3E). The VIP plot (Fig. 3F) describes tyramine, aspartate, asparagine, o-phosphocholine, glutamate, sn-Glycero-3-phosphocholine, proline, glycine, glutamine, lysine, valine and citrate as significant metabolites.

Pairwise comparison of sensory groups in Gouda type cheeses

The effect of variation on the ingredients of metabolite profiles and the concentration of specific metabolites in cheese will reveal valuable knowledge for product design and process development. In our study, we adopted a pairwise comparison in PLS-DA and OPLS-DA tools to obtain the VIP scores and a corresponding loading plot that outlines the significant metabolites in the discrimination of each sensory group. The pairwise comparison for the Gouda type cheese samples using PLS-DA models shows the clear separation of G1 against G2 and G3 samples (Fig. 4 A-C). Their corresponding VIP plots (Fig. 4 D-E), unveiled citrate as the most relevant metabolite involved in the discrimination. A partial separation between the G2 and G3 samples signifies the differences in metabolite profiles among them were not much better. The R2 and Q2 values were in the ranges of 0.430–0.748 and −0.284 to 0.447, respectively.

Fig. 4

PLS-DA score (A-C) and VIP (D-F) plots derived from the 1H NMR spectra of Gouda type cheese samples demonstrating the separation between the sensory groups of 1 and 2 (A, B), 1 and 3 (C, D) and 2 and 3 (E, F).

Similar, pairwise comparison on an OPLS-DA (Fig. 5 A-C) exhibits greater separation between the groups, than in PLS-DA. OPLS-DA provides a clear separation when variation between groups is greater than the variation within the group in the dataset. The model employs orthogonal signal correction (OSC) to filter the information that has no connection with the response matrix (Y) from the independent variable matrix (X), termed as “structured noise”. Thereby, OPLS-DA distinguishes the differences between groups with improved efficiency and analytical power (Worley and Powers, 2013). The OPLS-DA models for pairwise comparison of G1 with G2 and G3 samples illustrate the location of these groups on the opposite side of the Hotelling T2 ellipse, confirming higher variations in their metabolite profiles (Fig. 5 A-B). A partial overlapping prevails between G2 and G3 samples, as similar to PLS-DA (Fig. 5 C). The R2Y and Q2 values were in the ranges of 0.379–0.783 and −0.166 to 0.928), respectively. The loading plot of OPLS-DA models illustrates the promising metabolites or discriminant markers involved in the separation of the two groups. The p[1] and p(corr)[1] in the loading plot, are X- and Y-axes vectors of the predictive components, which define the significance and consistency of individual variables in X, respectively. The metabolites distant from the centre axis correspond to the greater deviation between the groups with more reliability and are relevant in the hunt for biomarkers (Worley and Powers, 2013). Our loading plots (Fig. 5 D-F) express lower levels of citrate, choline and O-phosphocholine in GI, higher levels of citrate in G3 and higher levels of lactate and citrate in G3 impel the separation of the sensory group among others. These data reveal the potential of these analytical and statistical approaches in characterizing Gouda type cheese samples; furthermore, they reveal citrate as the sole metabolite marker in discrimination of the sensory groups.

Fig. 5

OPLS-DA score plots (A-C) and corresponding loading plots (D-F) derived from the 1H NMR spectra of Gouda type cheese samples demonstrating the separation between the sensory groups of 1 and 2 (A, B), 1 and 3 (C, D) and 2 and 3 (E, F).

Pairwise comparison of sensory groups in other types cheeses

The PLS-DA models showed that samples of G1 are well distinguished from the other sensory groups in both G2 and G3 samples (Fig. 6 A-B). Furthermore, the separation was also observed between G2 and G3 samples (Fig. 6 C). These results confirmed the dissimilarity in the metabolite profiles between each sensory group in other types of cheeses. The VIP of the major metabolites owing to the separation of the sensory groups in the PLS-DA model is illustrated in Fig. 6 D-F. The R2 and Q2 values were in the ranges of 0.606–0.793 and −0.045 to 0.377, respectively. A further examination based on OPLS-DA score plots showed an improved separation with good reliability and predictive ability (Fig. 7). In both the PLS-DA and OPLS-DA models, the sensory groups were on the opposite side signifying a higher variation in their metabolite composition. The R2Y and Q2 values were in the ranges of 0.621–0.958 and 0.263–0.832, respectively. Compounds that contribute to the separation of sensory groups were extracted in the corresponding loading plots (Fig. 7 D-F). It is evident that the metabolite profiles are distinctly different in each sensory group. Higher levels of phospholipids (O-phosphocholine and sn-glycero-3-phosphocholine) and amino acids (glycine and asparagine) found in G I cheese samples greatly influenced their separation from G2 and G3. Likewise, higher levels of proline and glutamine in G2 samples and the metabolite tyramine in G3 cheese samples significantly influences the separation of these groups from others. These metabolites are also recognized as significant variables in the VIP plots of the corresponding pairwise comparison of sensory groups in PLS-DA models.

Fig. 6

PLS-DA score (A-C) and VIP (D-F) plots derived from the 1H NMR spectra of all cheese samples (except Gouda type) demonstrating the separation between the sensory groups of 1 and 2 (A, B), 1 and 3 (C, D) and 2 and 3 (E, F).

Fig. 7

OPLS-DA score plots (A-C) and corresponding loading plots (D-F) derived from the 1H NMR spectra of all cheese samples (except Gouda type) demonstrating the separation between the sensory groups of 1 and 2 (A, B), 1 and 3 (C, D) and 2 and 3 (E, F).

Heat map analysis of metabolite profiles

The quality and taste of the cheeses mainly depend upon the metabolite profiles and the precise balance of concentration of various metabolites. Identification and quantification of metabolites that significantly influence the sensory attributes would be much more informative. Hence, in our study we produced heatmaps to detect the unique profile of the cheeses and the metabolite variations among the sensory groups. The heatmap demonstrating the HR-MAS NMR metabolic profiles for each cheese sample is presented in Fig. 8A. The heatmap results are consistent with the results of the multivariate statistical analysis. The strong outlier (52-brown cheese) in the PCA (Fig. 2A), displays a quite distinct metabolite profile with higher levels of carbohydrates (lactose, galactose) and organic acids (citrate, pyruvate), absence of lactate and lower levels of amino acids. A similar NMR spectrum with higher levels of lactose and organic acids in probiotic cheeses was recorded at an early ripening stage (Rodrigues et al., 2011). Thus, the reduced lactose transformation and proteolytic activity observed in brown cheese might be associated with the early ripening conditions. Piras et al. (2013) demonstrated the rapid utilization of primary carbohydrates (lactose, galactose) by the starter and adjunct cultures in the earlier ripening stage (within 15 days). In accordance, lactose and galactose are detected only in brown cheese. Likewise, the weak outlier (36-Cheddar cheese) shows an elevated pyruvate level along with lower glycine, serine and proline levels (Fig. 2A). Pyruvate, the intermediate in lactose fermentation, was found to be prominent only in the outliers (brown and cheddar type cheeses). Furthermore, the presence of biogenic amine, tyramine in the Camembert (40) cheese alone was confirmed (Fig. 2).

Fig. 8

The HCA heatmap plots based on the 29 characteristic metabolites identified for all the cheese samples (A); clear metabolic differences between sensory groups in Gouda cheese type (B) and other cheese types (C). Colours from highest (red) to lowest (blue) represent metabolite expression values in different groups.

Metabolite profile discrimination of sensory groups

In this study, we also assessed the variation in the NMR metabolite levels of the cheese samples based on sensory groups. The heatmaps demonstrated clear discrimination between each sensory group in both the Gouda type (Fig. 8B) and other types of cheeses (Fig. 8C) indicating the metabolite composition of the cheese samples differs dramatically between the cheese samples in the different sensory groups. The cheese aroma characteristics are highly determined by the metabolites that yielded during proteolysis, metabolism of carbohydrates, and lipolysis. The detailed information about the starter cultures employed in the Korean dairy farms for the production of the cheese samples used in this study is summarised in Table 5. The starter culture bacteria employed in the fermentation process acts as a main supplier of the enzymes involved in the foresaid pathways. Predominantly, employment of Lactobacilli and Lactococci starters has been proven to hasten the cheese ripening process and enhancement of their organoleptic qualities (McSweeney and Sousa, 2000, Khattab et al., 2019). Proteolysis is the chief complex biochemical reaction that occurs during the cheese ripening process, responsible for the metabolism of milk proteins into small peptides and free amino acids by the proteolytic enzymes of starter cultures. The amino acids content acts as an indicator of the proteolytic activity and for higher contribution in the background texture and flavour of several cheese varieties (Ruyssen et al., 2013, Khattab et al., 2019). In our study, the proportion of free amino acids was significantly higher in all types of cheese samples in G1 comparatively with G2 or G3. These results were in agreement with earlier studies (Scano et al., 2019b, Rodrigues et al., 2011) suggesting that cheese samples in G1 underwent a more extensive proteolysis process. A higher degree of proteolytic activity during the initial ripening period resulted in increased peptide concentration that favours the microbial activity of aminopeptidases/exopeptidases resulting in the evolution of free amino acids during the second period of ripening (30–60 days) (Rodrigues et al., 2011). Moreover, the levels of umami tasting amino acids such as glutamate, aspartate and asparagine that contribute to the taste of cheese were also considerably higher in G1 than the other two groups. This could mean that the proceeded ripening periods increased the amino acids content along with subsequent flavour development by amino acids themselves and as precursors of volatile compounds that lead to the strong taste and aroma while reducing the bitterness in G1 cheese samples.

Table 5

Starter strains employed for the production of cheese samples used in this study.

Taxonomy	Culture (Company)	Sample nos.	Cheese type
Lactococcus lactis subsp. lactis,L. lactis subsp. cremoris,L. lactis subsp. diacetylactis	MM100(Danisco)	1,9,14,18,25,28,38,40,44	Gouda
Lactococcus lactis subsp. cremoris,L. lactis subsp. diacetylactis,L. lactis subsp. lactis,Leuconostoc mensenteroides subsp.cremoris	CHN-11(Chr.Hansen)	3,4,5,6,7,8,11,12,19,21,29,31,36, 37,41,42,43,45, 49,51	Gouda
		13,24,27,39,47	Camembert
		48	Brie
	MW 046N(SACCO)	16,30,34,50	Gouda
	Provat 322(Danisco)	10,35	Cheddar
		32	Quark
		46	Gouda
Lactobacillus delbrueckii subsp. bulgaricus,Streptococcus thermophilus	Y 083F(SACCO)	2	Berg
Lactobacillus bulgaricus,Steptococcus thermophilus	TCC-3(Chr.Hansen)	17
		26	Appenzeller
Streptococcus thermophilus,Lactobacillus delbrueckii subsp. lactis,Lactobacillus helveticus	SuCasu(Danisco)	22	Berg
		33	Appenzeller
Lactobacillus bulgaricus,Steptococcus thermophilus,Lactobacillus helveticus	TCC-3,LB-B02(Chr.Hansen)	15	Frill
Lactobacillus bulgaricus,Steptococcus thermophilus,Lactococcus lactis subsp. cremoris,L. lactis subsp. diacetylactis,L. lactis subsp. lactis,Leuconostoc mensenteroides subsp. cremoris	TCC-3,CHN-11(Chr.Hansen)	20	Reblochon
Streptococcus thermophilus,Lactobacillus delbrueckii subsp. lactis,L. helveticus,Propionibacterium freudenreichii subsp.shermanii	SuCasu,Choozit Eyes 2(Danisco)	23	Emmental

In conflict, the relative amount of citrate was lower in G1 in comparison with the other two groups. These results could be correlated to the metabolism of citrate into pyruvate and acetate by the citrate positive (Cit+) strains. During the cheese-making process, most of the soluble citrate (90%) in fresh milk is removed as whey, while the retained citrate is metabolised as diverse flavour compounds (acetate, 2,3-diacetyl and acetoin). Since the starters (Cit+) cannot use citrate as an energy source, co-metabolism of citrate occurs in the presence of fermentable carbohydrates (Fox et al., 2017). Furthermore, the CO2 produced during the citrate metabolism influence the texture imparting the “eye” structure in some Gouda cheeses (McSweeney and Sousa, 2000). Studies on the association between the bovine metabolites and rennet coagulation properties of milk showed an increasing trend of the citrate concentration in samples with poor coagulation properties (Sundekilde et al., 2011).

Similar to citrate, the most dominant metabolite lactate was in higher levels in G3 than the other two groups. The rapid transformation of lactose to lactate during the initial ripening process of cheese is crucial in the production of all cheese varieties. Lactate certainly promotes the flavour in acid-curd and ripened cheese varieties via transforming glycolysis or phosphoketolase pathways, depending upon the utilization of starter cultures (Jo et al., 2018). Other metabolites such as glycerol, choline and phospholipids (O-phosphocholine and sn-Glycero-3-phosphocholine) were relatively in a lower level in sensory G1 of Gouda type cheeses. In contrast, these metabolites were relatively high in GI cheese samples from other types of cheeses. The combination of metabolite profiling and multivariate analysis delivered significant information on the composition and concentration of metabolites in sensory group discrimination of cheeses.

Conclusion

The current study demonstrates 1H HRMAS-NMR based metabolic fingerprinting as a valuable tool in distinguishing cheese samples. The NMR data combined with chemometric analysis considerably improves the discrimination of groups based on sensory qualities and reveals potential biomarkers. Interestingly, our studies showed the NMR based metabolite profiles of the domestic Korean cheeses are in agreement with their sensory characteristics. In contrary to earlier studies, the location of cheese samples in the multivariate models was based on NMR data set, rather than on their geographical origin. Amino acids, lactate, citrate and phospholipids were the most relevant variables/markers related to their sensory qualities. The metabolite variations in the cheeses are likely to be influenced by the degree of proteolytic and metabolic activities that occur during cheese ripening. Further, profound studies on investigating the consistency of metabolic fingerprints could lead to the establishment of discrimination biomarkers for sensory determination in the cheese industry. Reliable discrimination of biomarkers associated with sensory traits will benefit the cheese industries to improve the taste and quality.

Declaration of Competing Interest

The authors declare that they have no known competing financial interests or personal relationships that could have appeared to influence the work reported in this paper.