Project Overview:

Objective: The goal is to use mathematical and statistical techniques to analyze milk yield in dairy cows, optimize feed intake, and predict trends in cow health using synthetic or open-source data. The project will demonstrate skills in linear algebra, statistics, and mathematical modeling using R.

# Load necessary libraries
library(ggplot2)
library(dplyr)
library(tidyr)
library(MASS)  # For synthetic data generation
library(stats)  # For PCA and ANOVA
library(caret)  # For Linear Regression
library(corrplot)  # For correlation matrix visualization
corrplot 0.95 loaded
# 1.1 Load open-source data (if available) from Kaggle or FAO
# Here we simulate the data as an example using synthetic data

# Synthetic data generation for cow feed intake and milk yield
set.seed(123)

# Create synthetic data for 100 cows over 365 days
n_cows <- 100
n_days <- 365
cow_ids <- rep(1:n_cows, each = n_days)

# Generate synthetic feed intake data (in kg) and milk yield (liters)
feed_protein <- rnorm(n_cows * n_days, mean = 16, sd = 2)
feed_fat <- rnorm(n_cows * n_days, mean = 4, sd = 1)
feed_energy <- rnorm(n_cows * n_days, mean = 2800, sd = 500)  # in kcal

milk_yield <- 25 + 0.5 * feed_protein + 0.3 * feed_fat + rnorm(n_cows * n_days, mean = 0, sd = 5)  # Milk yield in liters

# Combine into a data frame
synthetic_data <- data.frame(
  cow_id = cow_ids,
  feed_protein = feed_protein,
  feed_fat = feed_fat,
  feed_energy = feed_energy,
  milk_yield = milk_yield
)

# View the first few rows
head(synthetic_data)
ABCDEFGHIJ0123456789
 
 
cow_id
<int>
feed_protein
<dbl>
feed_fat
<dbl>
feed_energy
<dbl>
milk_yield
<dbl>
1114.879052.5108112848.70228.20588
2115.539653.8954534234.10636.41569
3119.117424.1255292889.79731.95770
4116.141022.9517162543.28533.59276
5116.258582.9669572142.22438.55488
6119.430134.3595331981.63133.60346
6 rows

# Summary of the synthetic data
summary(synthetic_data)
     cow_id        feed_protein       feed_fat      
 Min.   :  1.00   Min.   : 8.309   Min.   :-0.1291  
 1st Qu.: 25.75   1st Qu.:14.646   1st Qu.: 3.3328  
 Median : 50.50   Median :15.978   Median : 4.0118  
 Mean   : 50.50   Mean   :15.988   Mean   : 4.0069  
 3rd Qu.: 75.25   3rd Qu.:17.348   3rd Qu.: 4.6770  
 Max.   :100.00   Max.   :24.646   Max.   : 7.8689  
  feed_energy     milk_yield   
 Min.   : 609   Min.   :13.36  
 1st Qu.:2469   1st Qu.:30.76  
 Median :2803   Median :34.17  
 Mean   :2803   Mean   :34.20  
 3rd Qu.:3138   3rd Qu.:37.63  
 Max.   :4903   Max.   :54.57

Step 2: Exploratory Data Analysis (EDA)

# 2.1 Summary statistics
summary(synthetic_data)
     cow_id        feed_protein       feed_fat        feed_energy     milk_yield   
 Min.   :  1.00   Min.   : 8.309   Min.   :-0.1291   Min.   : 609   Min.   :13.36  
 1st Qu.: 25.75   1st Qu.:14.646   1st Qu.: 3.3328   1st Qu.:2469   1st Qu.:30.76  
 Median : 50.50   Median :15.978   Median : 4.0118   Median :2803   Median :34.17  
 Mean   : 50.50   Mean   :15.988   Mean   : 4.0069   Mean   :2803   Mean   :34.20  
 3rd Qu.: 75.25   3rd Qu.:17.348   3rd Qu.: 4.6770   3rd Qu.:3138   3rd Qu.:37.63  
 Max.   :100.00   Max.   :24.646   Max.   : 7.8689   Max.   :4903   Max.   :54.57  
# 2.2 Correlation matrix
cor_matrix <- cor(synthetic_data[, c("feed_protein", "feed_fat", "feed_energy", "milk_yield")])
print(cor_matrix)
             feed_protein     feed_fat feed_energy  milk_yield
feed_protein  1.000000000 -0.002192113 0.002754189 0.195083280
feed_fat     -0.002192113  1.000000000 0.001116323 0.053837849
feed_energy   0.002754189  0.001116323 1.000000000 0.004458438
milk_yield    0.195083280  0.053837849 0.004458438 1.000000000
# Visualize the correlation matrix
corrplot(cor_matrix, method = "color", type = "upper")


# 2.3 Visualizations
# Plot milk yield vs feed protein intake
ggplot(synthetic_data, aes(x = feed_protein, y = milk_yield)) +
  geom_point(alpha = 0.5) +
  geom_smooth(method = "lm", color = "blue") +
  labs(title = "Milk Yield vs Protein Intake", x = "Protein Intake (kg)", y = "Milk Yield (liters)")


# Plot milk yield vs feed fat intake
ggplot(synthetic_data, aes(x = feed_fat, y = milk_yield)) +
  geom_point(alpha = 0.5) +
  geom_smooth(method = "lm", color = "red") +
  labs(title = "Milk Yield vs Fat Intake", x = "Fat Intake (kg)", y = "Milk Yield (liters)")

NA
NA

Step 3: Linear Regression for Milk Yield Prediction

We will use linear regression to predict milk yield based on the feed intake of protein, fat, and energy.

# 3.1 Linear regression: Milk yield as a function of feed intake
lm_model <- lm(milk_yield ~ feed_protein + feed_fat + feed_energy, data = synthetic_data)
summary(lm_model)

Call:
lm(formula = milk_yield ~ feed_protein + feed_fat + feed_energy, 
    data = synthetic_data)

Residuals:
     Min       1Q   Median       3Q      Max 
-20.9205  -3.3634  -0.0308   3.3950  20.0594 

Coefficients:
              Estimate Std. Error t value Pr(>|t|)    
(Intercept)  2.499e+01  2.778e-01  89.961   <2e-16 ***
feed_protein 4.999e-01  1.313e-02  38.078   <2e-16 ***
feed_fat     2.767e-01  2.614e-02  10.585   <2e-16 ***
feed_energy  3.943e-05  5.237e-05   0.753    0.451    
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Residual standard error: 5.008 on 36496 degrees of freedom
Multiple R-squared:  0.04102,   Adjusted R-squared:  0.04094 
F-statistic: 520.3 on 3 and 36496 DF,  p-value: < 2.2e-16
# Visualize regression diagnostics
par(mfrow = c(2, 2))
plot(lm_model)


# 3.2 Predictions
predicted_milk_yield <- predict(lm_model, synthetic_data)

# Add the predicted values to the dataset
synthetic_data$predicted_milk_yield <- predicted_milk_yield

# 3.3 Visualize the actual vs predicted milk yield
ggplot(synthetic_data, aes(x = milk_yield, y = predicted_milk_yield)) +
  geom_point(alpha = 0.5, color = "green") +
  geom_abline(slope = 1, intercept = 0, color = "red") +
  labs(title = "Actual vs Predicted Milk Yield", x = "Actual Milk Yield (liters)", y = "Predicted Milk Yield (liters)")

Step 4: Principal Component Analysis (PCA)

Perform PCA on the feed intake data to identify the most important features contributing to milk yield.

# 4.1 Perform PCA on feed intake data
pca_data <- synthetic_data[, c("feed_protein", "feed_fat", "feed_energy")]
pca_result <- prcomp(pca_data, scale. = TRUE)

# 4.2 Summary of PCA results
summary(pca_result)
Importance of components:
                          PC1    PC2    PC3
Standard deviation     1.0015 1.0005 0.9979
Proportion of Variance 0.3343 0.3337 0.3320
Cumulative Proportion  0.3343 0.6680 1.0000
# 4.3 Visualize PCA (Scree plot)
screeplot(pca_result, type = "lines", main = "Scree Plot")


# 4.4 Biplot of PCA
biplot(pca_result)

Step 5: Linear Algebra in Nutrient Calculation

We will use matrix operations to calculate the total nutrient intake for a group of cows.

# 5.1 Matrix operations for total nutrient intake
nutrient_matrix <- as.matrix(synthetic_data[, c("feed_protein", "feed_fat", "feed_energy")])
nutrient_totals <- colSums(nutrient_matrix)

# 5.2 Display the total nutrient intake
print(nutrient_totals)
feed_protein     feed_fat  feed_energy 
    583553.2     146252.0  102321116.6

Step 6: Statistical Testing (ANOVA)

We will test whether there are significant differences in milk yield between different levels of protein intake using ANOVA.

# 6.1 Categorize protein intake into low, medium, high categories
synthetic_data$protein_level <- cut(synthetic_data$feed_protein, breaks = 3, labels = c("Low", "Medium", "High"))

# 6.2 Perform ANOVA on milk yield by protein level
anova_result <- aov(milk_yield ~ protein_level, data = synthetic_data)
summary(anova_result)
                 Df Sum Sq Mean Sq F value Pr(>F)    
protein_level     2  20932   10466   409.2 <2e-16 ***
Residuals     36497 933541      26                   
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
# 6.3 Visualize milk yield by protein level
ggplot(synthetic_data, aes(x = protein_level, y = milk_yield, fill = protein_level)) +
  geom_boxplot() +
  labs(title = "Milk Yield by Protein Level", x = "Protein Level", y = "Milk Yield (liters)")

LS0tCnRpdGxlOiAiTGl2ZXN0b2NrIEZlZWQgYW5kIE1pbGsgWWllbGQgQW5hbHlzaXMgd2l0aCBSIgpvdXRwdXQ6CiAgaHRtbF9ub3RlYm9vazogZGVmYXVsdAogIHBkZl9kb2N1bWVudDogZGVmYXVsdAotLS0KCiMgUHJvamVjdCBPdmVydmlldzoKT2JqZWN0aXZlOiBUaGUgZ29hbCBpcyB0byB1c2UgbWF0aGVtYXRpY2FsIGFuZCBzdGF0aXN0aWNhbCB0ZWNobmlxdWVzIHRvIGFuYWx5emUgbWlsayB5aWVsZCBpbiBkYWlyeSBjb3dzLCBvcHRpbWl6ZSBmZWVkIGludGFrZSwgYW5kIHByZWRpY3QgdHJlbmRzIGluIGNvdyBoZWFsdGggdXNpbmcgc3ludGhldGljIG9yIG9wZW4tc291cmNlIGRhdGEuIFRoZSBwcm9qZWN0IHdpbGwgZGVtb25zdHJhdGUgc2tpbGxzIGluIGxpbmVhciBhbGdlYnJhLCBzdGF0aXN0aWNzLCBhbmQgbWF0aGVtYXRpY2FsIG1vZGVsaW5nIHVzaW5nIFIuCgpgYGB7cn0KIyBMb2FkIG5lY2Vzc2FyeSBsaWJyYXJpZXMKbGlicmFyeShnZ3Bsb3QyKQpsaWJyYXJ5KGRwbHlyKQpsaWJyYXJ5KHRpZHlyKQpsaWJyYXJ5KE1BU1MpICAjIEZvciBzeW50aGV0aWMgZGF0YSBnZW5lcmF0aW9uCmxpYnJhcnkoc3RhdHMpICAjIEZvciBQQ0EgYW5kIEFOT1ZBCmxpYnJhcnkoY2FyZXQpICAjIEZvciBMaW5lYXIgUmVncmVzc2lvbgpsaWJyYXJ5KGNvcnJwbG90KSAgIyBGb3IgY29ycmVsYXRpb24gbWF0cml4IHZpc3VhbGl6YXRpb24KCiMgMS4xIExvYWQgb3Blbi1zb3VyY2UgZGF0YSAoaWYgYXZhaWxhYmxlKSBmcm9tIEthZ2dsZSBvciBGQU8KIyBIZXJlIHdlIHNpbXVsYXRlIHRoZSBkYXRhIGFzIGFuIGV4YW1wbGUgdXNpbmcgc3ludGhldGljIGRhdGEKCiMgU3ludGhldGljIGRhdGEgZ2VuZXJhdGlvbiBmb3IgY293IGZlZWQgaW50YWtlIGFuZCBtaWxrIHlpZWxkCnNldC5zZWVkKDEyMykKCiMgQ3JlYXRlIHN5bnRoZXRpYyBkYXRhIGZvciAxMDAgY293cyBvdmVyIDM2NSBkYXlzCm5fY293cyA8LSAxMDAKbl9kYXlzIDwtIDM2NQpjb3dfaWRzIDwtIHJlcCgxOm5fY293cywgZWFjaCA9IG5fZGF5cykKCiMgR2VuZXJhdGUgc3ludGhldGljIGZlZWQgaW50YWtlIGRhdGEgKGluIGtnKSBhbmQgbWlsayB5aWVsZCAobGl0ZXJzKQpmZWVkX3Byb3RlaW4gPC0gcm5vcm0obl9jb3dzICogbl9kYXlzLCBtZWFuID0gMTYsIHNkID0gMikKZmVlZF9mYXQgPC0gcm5vcm0obl9jb3dzICogbl9kYXlzLCBtZWFuID0gNCwgc2QgPSAxKQpmZWVkX2VuZXJneSA8LSBybm9ybShuX2Nvd3MgKiBuX2RheXMsIG1lYW4gPSAyODAwLCBzZCA9IDUwMCkgICMgaW4ga2NhbAoKbWlsa195aWVsZCA8LSAyNSArIDAuNSAqIGZlZWRfcHJvdGVpbiArIDAuMyAqIGZlZWRfZmF0ICsgcm5vcm0obl9jb3dzICogbl9kYXlzLCBtZWFuID0gMCwgc2QgPSA1KSAgIyBNaWxrIHlpZWxkIGluIGxpdGVycwoKIyBDb21iaW5lIGludG8gYSBkYXRhIGZyYW1lCnN5bnRoZXRpY19kYXRhIDwtIGRhdGEuZnJhbWUoCiAgY293X2lkID0gY293X2lkcywKICBmZWVkX3Byb3RlaW4gPSBmZWVkX3Byb3RlaW4sCiAgZmVlZF9mYXQgPSBmZWVkX2ZhdCwKICBmZWVkX2VuZXJneSA9IGZlZWRfZW5lcmd5LAogIG1pbGtfeWllbGQgPSBtaWxrX3lpZWxkCikKCiMgVmlldyB0aGUgZmlyc3QgZmV3IHJvd3MKaGVhZChzeW50aGV0aWNfZGF0YSkKCiMgU3VtbWFyeSBvZiB0aGUgc3ludGhldGljIGRhdGEKc3VtbWFyeShzeW50aGV0aWNfZGF0YSkKCmBgYAoKCiMjIFN0ZXAgMjogRXhwbG9yYXRvcnkgRGF0YSBBbmFseXNpcyAoRURBKQoKYGBge3J9CiMgMi4xIFN1bW1hcnkgc3RhdGlzdGljcwpzdW1tYXJ5KHN5bnRoZXRpY19kYXRhKQoKIyAyLjIgQ29ycmVsYXRpb24gbWF0cml4CmNvcl9tYXRyaXggPC0gY29yKHN5bnRoZXRpY19kYXRhWywgYygiZmVlZF9wcm90ZWluIiwgImZlZWRfZmF0IiwgImZlZWRfZW5lcmd5IiwgIm1pbGtfeWllbGQiKV0pCnByaW50KGNvcl9tYXRyaXgpCgojIFZpc3VhbGl6ZSB0aGUgY29ycmVsYXRpb24gbWF0cml4CmNvcnJwbG90KGNvcl9tYXRyaXgsIG1ldGhvZCA9ICJjb2xvciIsIHR5cGUgPSAidXBwZXIiKQoKIyAyLjMgVmlzdWFsaXphdGlvbnMKIyBQbG90IG1pbGsgeWllbGQgdnMgZmVlZCBwcm90ZWluIGludGFrZQpnZ3Bsb3Qoc3ludGhldGljX2RhdGEsIGFlcyh4ID0gZmVlZF9wcm90ZWluLCB5ID0gbWlsa195aWVsZCkpICsKICBnZW9tX3BvaW50KGFscGhhID0gMC41KSArCiAgZ2VvbV9zbW9vdGgobWV0aG9kID0gImxtIiwgY29sb3IgPSAiYmx1ZSIpICsKICBsYWJzKHRpdGxlID0gIk1pbGsgWWllbGQgdnMgUHJvdGVpbiBJbnRha2UiLCB4ID0gIlByb3RlaW4gSW50YWtlIChrZykiLCB5ID0gIk1pbGsgWWllbGQgKGxpdGVycykiKQoKIyBQbG90IG1pbGsgeWllbGQgdnMgZmVlZCBmYXQgaW50YWtlCmdncGxvdChzeW50aGV0aWNfZGF0YSwgYWVzKHggPSBmZWVkX2ZhdCwgeSA9IG1pbGtfeWllbGQpKSArCiAgZ2VvbV9wb2ludChhbHBoYSA9IDAuNSkgKwogIGdlb21fc21vb3RoKG1ldGhvZCA9ICJsbSIsIGNvbG9yID0gInJlZCIpICsKICBsYWJzKHRpdGxlID0gIk1pbGsgWWllbGQgdnMgRmF0IEludGFrZSIsIHggPSAiRmF0IEludGFrZSAoa2cpIiwgeSA9ICJNaWxrIFlpZWxkIChsaXRlcnMpIikKCgpgYGAKCiMjIFN0ZXAgMzogTGluZWFyIFJlZ3Jlc3Npb24gZm9yIE1pbGsgWWllbGQgUHJlZGljdGlvbgpXZSB3aWxsIHVzZSBsaW5lYXIgcmVncmVzc2lvbiB0byBwcmVkaWN0IG1pbGsgeWllbGQgYmFzZWQgb24gdGhlIGZlZWQgaW50YWtlIG9mIHByb3RlaW4sIGZhdCwgYW5kIGVuZXJneS4KCmBgYHtyfQojIDMuMSBMaW5lYXIgcmVncmVzc2lvbjogTWlsayB5aWVsZCBhcyBhIGZ1bmN0aW9uIG9mIGZlZWQgaW50YWtlCmxtX21vZGVsIDwtIGxtKG1pbGtfeWllbGQgfiBmZWVkX3Byb3RlaW4gKyBmZWVkX2ZhdCArIGZlZWRfZW5lcmd5LCBkYXRhID0gc3ludGhldGljX2RhdGEpCnN1bW1hcnkobG1fbW9kZWwpCgojIFZpc3VhbGl6ZSByZWdyZXNzaW9uIGRpYWdub3N0aWNzCnBhcihtZnJvdyA9IGMoMiwgMikpCnBsb3QobG1fbW9kZWwpCgojIDMuMiBQcmVkaWN0aW9ucwpwcmVkaWN0ZWRfbWlsa195aWVsZCA8LSBwcmVkaWN0KGxtX21vZGVsLCBzeW50aGV0aWNfZGF0YSkKCiMgQWRkIHRoZSBwcmVkaWN0ZWQgdmFsdWVzIHRvIHRoZSBkYXRhc2V0CnN5bnRoZXRpY19kYXRhJHByZWRpY3RlZF9taWxrX3lpZWxkIDwtIHByZWRpY3RlZF9taWxrX3lpZWxkCgojIDMuMyBWaXN1YWxpemUgdGhlIGFjdHVhbCB2cyBwcmVkaWN0ZWQgbWlsayB5aWVsZApnZ3Bsb3Qoc3ludGhldGljX2RhdGEsIGFlcyh4ID0gbWlsa195aWVsZCwgeSA9IHByZWRpY3RlZF9taWxrX3lpZWxkKSkgKwogIGdlb21fcG9pbnQoYWxwaGEgPSAwLjUsIGNvbG9yID0gImdyZWVuIikgKwogIGdlb21fYWJsaW5lKHNsb3BlID0gMSwgaW50ZXJjZXB0ID0gMCwgY29sb3IgPSAicmVkIikgKwogIGxhYnModGl0bGUgPSAiQWN0dWFsIHZzIFByZWRpY3RlZCBNaWxrIFlpZWxkIiwgeCA9ICJBY3R1YWwgTWlsayBZaWVsZCAobGl0ZXJzKSIsIHkgPSAiUHJlZGljdGVkIE1pbGsgWWllbGQgKGxpdGVycykiKQoKYGBgCgojIyBTdGVwIDQ6IFByaW5jaXBhbCBDb21wb25lbnQgQW5hbHlzaXMgKFBDQSkKUGVyZm9ybSBQQ0Egb24gdGhlIGZlZWQgaW50YWtlIGRhdGEgdG8gaWRlbnRpZnkgdGhlIG1vc3QgaW1wb3J0YW50IGZlYXR1cmVzIGNvbnRyaWJ1dGluZyB0byBtaWxrIHlpZWxkLgoKYGBge3J9CiMgNC4xIFBlcmZvcm0gUENBIG9uIGZlZWQgaW50YWtlIGRhdGEKcGNhX2RhdGEgPC0gc3ludGhldGljX2RhdGFbLCBjKCJmZWVkX3Byb3RlaW4iLCAiZmVlZF9mYXQiLCAiZmVlZF9lbmVyZ3kiKV0KcGNhX3Jlc3VsdCA8LSBwcmNvbXAocGNhX2RhdGEsIHNjYWxlLiA9IFRSVUUpCgojIDQuMiBTdW1tYXJ5IG9mIFBDQSByZXN1bHRzCnN1bW1hcnkocGNhX3Jlc3VsdCkKCiMgNC4zIFZpc3VhbGl6ZSBQQ0EgKFNjcmVlIHBsb3QpCnNjcmVlcGxvdChwY2FfcmVzdWx0LCB0eXBlID0gImxpbmVzIiwgbWFpbiA9ICJTY3JlZSBQbG90IikKCiMgNC40IEJpcGxvdCBvZiBQQ0EKYmlwbG90KHBjYV9yZXN1bHQpCgpgYGAKCiMjIFN0ZXAgNTogTGluZWFyIEFsZ2VicmEgaW4gTnV0cmllbnQgQ2FsY3VsYXRpb24KV2Ugd2lsbCB1c2UgbWF0cml4IG9wZXJhdGlvbnMgdG8gY2FsY3VsYXRlIHRoZSB0b3RhbCBudXRyaWVudCBpbnRha2UgZm9yIGEgZ3JvdXAgb2YgY293cy4KCmBgYHtyfQojIDUuMSBNYXRyaXggb3BlcmF0aW9ucyBmb3IgdG90YWwgbnV0cmllbnQgaW50YWtlCm51dHJpZW50X21hdHJpeCA8LSBhcy5tYXRyaXgoc3ludGhldGljX2RhdGFbLCBjKCJmZWVkX3Byb3RlaW4iLCAiZmVlZF9mYXQiLCAiZmVlZF9lbmVyZ3kiKV0pCm51dHJpZW50X3RvdGFscyA8LSBjb2xTdW1zKG51dHJpZW50X21hdHJpeCkKCiMgNS4yIERpc3BsYXkgdGhlIHRvdGFsIG51dHJpZW50IGludGFrZQpwcmludChudXRyaWVudF90b3RhbHMpCgpgYGAKCiMjIFN0ZXAgNjogU3RhdGlzdGljYWwgVGVzdGluZyAoQU5PVkEpCldlIHdpbGwgdGVzdCB3aGV0aGVyIHRoZXJlIGFyZSBzaWduaWZpY2FudCBkaWZmZXJlbmNlcyBpbiBtaWxrIHlpZWxkIGJldHdlZW4gZGlmZmVyZW50IGxldmVscyBvZiBwcm90ZWluIGludGFrZSB1c2luZyBBTk9WQS4KCmBgYHtyfQojIDYuMSBDYXRlZ29yaXplIHByb3RlaW4gaW50YWtlIGludG8gbG93LCBtZWRpdW0sIGhpZ2ggY2F0ZWdvcmllcwpzeW50aGV0aWNfZGF0YSRwcm90ZWluX2xldmVsIDwtIGN1dChzeW50aGV0aWNfZGF0YSRmZWVkX3Byb3RlaW4sIGJyZWFrcyA9IDMsIGxhYmVscyA9IGMoIkxvdyIsICJNZWRpdW0iLCAiSGlnaCIpKQoKIyA2LjIgUGVyZm9ybSBBTk9WQSBvbiBtaWxrIHlpZWxkIGJ5IHByb3RlaW4gbGV2ZWwKYW5vdmFfcmVzdWx0IDwtIGFvdihtaWxrX3lpZWxkIH4gcHJvdGVpbl9sZXZlbCwgZGF0YSA9IHN5bnRoZXRpY19kYXRhKQpzdW1tYXJ5KGFub3ZhX3Jlc3VsdCkKCiMgNi4zIFZpc3VhbGl6ZSBtaWxrIHlpZWxkIGJ5IHByb3RlaW4gbGV2ZWwKZ2dwbG90KHN5bnRoZXRpY19kYXRhLCBhZXMoeCA9IHByb3RlaW5fbGV2ZWwsIHkgPSBtaWxrX3lpZWxkLCBmaWxsID0gcHJvdGVpbl9sZXZlbCkpICsKICBnZW9tX2JveHBsb3QoKSArCiAgbGFicyh0aXRsZSA9ICJNaWxrIFlpZWxkIGJ5IFByb3RlaW4gTGV2ZWwiLCB4ID0gIlByb3RlaW4gTGV2ZWwiLCB5ID0gIk1pbGsgWWllbGQgKGxpdGVycykiKQoKYGBgCgo=