In ScottHMcKean/duvernay_geomech_model: A Geostatistical Model of Geomechanical Properties in the Duvernay Play

library(duvernaygeomechmodel)

Load datasets

mldf_wres = read_csv("../data/merged_tests_w_residual.csv") %>%
  mutate(brittleness = (peak_stress_dev - residual_stress_dev)/peak_stress_dev) %>%
  mutate(orientation = as_factor(orientation)) %>%
  mutate(group = as_factor(group)) %>%
  rename(
    "Confining Pressure (MPa)" = confining_pressure,
    "Peak Stress (MPa)" = peak_stress_dev, 
    "Residual Stress (MPa)" = residual_stress_dev, 
    "Young's Modulus (GPa)" = youngs_modulus,
    "Poisson's Ratio" = poisson_ratio, 
    "Bulk Density (g/cm3)" = bulk_density, 
    "Brittleness (MPa)" = brittleness, 
    "Compressional Sonic (m/s)" = pred_vp, 
    "Shear Sonic (m/s)" = pred_vs,
    "Orientation" = orientation,
    "Group" = group)

mldf_nores = read_csv("../data/merged_tests_no_residual.csv") %>%
  mutate(orientation = as_factor(orientation)) %>%
  mutate(group = as_factor(group)) %>%
  rename(
    "Confining Pressure (MPa)" = confining_pressure,
    "Peak Stress (MPa)" = peak_stress_dev, 
    "Young's Modulus (GPa)" = youngs_modulus,
    "Poisson's Ratio" = poisson_ratio, 
    "Bulk Density (g/cm3)" = bulk_density,
    "Compressional Sonic (m/s)" = pred_vp, 
    "Shear Sonic (m/s)" = pred_vs,
    "Orientation" = orientation,
    "Group" = group)

mldf_wres_num <- mldf_wres %>% 
  select("Confining Pressure (MPa)", "Peak Stress (MPa)", 
    "Residual Stress (MPa)", "Young's Modulus (GPa)", "Poisson's Ratio", 
    "Bulk Density (g/cm3)", "Compressional Sonic (m/s)", "Shear Sonic (m/s)", 
    "Brittleness (MPa)")

mldf_nores_num <- mldf_nores %>% 
  select("Confining Pressure (MPa)", "Peak Stress (MPa)", 
    "Young's Modulus (GPa)", "Poisson's Ratio", "Bulk Density (g/cm3)", 
    "Compressional Sonic (m/s)", "Shear Sonic (m/s)")

Correlation Plot

cor_nores <- ggcorrplot(cor(mldf_nores_num),lab = TRUE, type = 'lower') +
  theme(legend.position = c(0.2,0.75)) +
  ggtitle("No Residual Stress (n = 410)")

cor_wres <- ggcorrplot(cor(mldf_wres_num),lab = TRUE, type = 'lower') +
  theme(legend.position = c(0.2,0.75)) +
  ggtitle("With Residual Stress (n = 304)")

grobs = arrangeGrob(grobs=list(cor_nores, cor_wres), widths = unit(c(6,6),"in"), heights = unit(c(6),"in"))
ggsave("../output/corrplot.jpg", grobs, width = 12, height = 6)
ggsave("../output/corrplot.pdf", grobs, width = 12, height = 6, dpi=600)

QQPlot

mldf_no_res_gather <- mldf_nores_num %>%
  gather()

mldf_wres_gather <- mldf_wres_num %>%
  gather()

ggplot(mldf_no_res_gather, aes(sample = value)) +
  geom_qq(size=1) +
  geom_qq_line() +
  geom_qq(data = mldf_wres_gather, aes(sample = value), color='red', size=1) +
  geom_qq_line(data = mldf_wres_gather, aes(sample = value), color='red') +
  facet_wrap(. ~ key, scales = 'free_y') +
  ylab('Sample Value') +
  xlab('Theoretical Value') +
  ggsave('../output/qq_plot_facet.pdf', height = 6, width = 12) +
  ggsave('../output/qq_plot_facet.pdf', height = 6, width = 12, dpi=600)

Box Plots

ggplot(mldf_no_res_gather, aes(y = value)) +
  geom_violin(aes(x='No Residual'), color='black', fill='grey') +
  geom_violin(data = mldf_wres_gather, aes(y = value, x='With Residual'), color='red', fill='grey') +
  geom_boxplot(aes(x='No Residual'), color='black', alpha=0.5) +
  geom_boxplot(data = mldf_wres_gather, aes(y = value, x='With Residual'), color='red', alpha=0.5) +
  facet_wrap(. ~ key, scales = 'free_y')

Filter-based feature importance using MLR and information gain - inconclusive.

# With Residuals
wres_feats = c("Confining Pressure (MPa)", "Peak Stress (MPa)", 
    "Residual Stress (MPa)", "Young's Modulus (GPa)", "Poisson's Ratio", 
    "Bulk Density (g/cm3)", "Compressional Sonic (m/s)", "Shear Sonic (m/s)", 
    "Brittleness (MPa)","Orientation","Group")

mldf_wres_clean <- mldf_wres %>% 
  select(all_of(wres_feats)) %>% 
  janitor::clean_names("upper_camel")

mldf_wres_ym_norm = normalizeFeatures(
  mldf_wres_clean %>% select(-PoissonsRatio, -BrittlenessMPa), 
  target="YoungsModulusGPa")

regr_wres_fv_ym = generateFilterValuesData(
  makeRegrTask(data = mldf_wres_ym_norm, target = "YoungsModulusGPa"), 
  method = c("FSelector_information.gain")
  )

ym_fv_wres_plot <- plotFilterValues(regr_wres_fv_ym) +
  ggtitle('YM w/ Residual Stress') +
  ylim(0,1) +
  coord_flip()

mldf_wres_pr_norm = normalizeFeatures(
  mldf_wres_clean %>% select(-YoungsModulusGPa, -BrittlenessMPa),
  target="PoissonsRatio")

regr_fv_wres_pr = generateFilterValuesData(
  makeRegrTask(data = mldf_wres_pr_norm, target = "PoissonsRatio"), 
  method = c("FSelector_information.gain")
  )

pr_fv_wres_plot <- plotFilterValues(regr_fv_wres_pr) +
  ggtitle('PR w/ Residual Stress') +
  ylim(0,1) +
  coord_flip()

mldf_wres_br_norm = normalizeFeatures(
  mldf_wres_clean %>% select(-YoungsModulusGPa, -PoissonsRatio), 
  target="BrittlenessMPa")

regr_fv_br = generateFilterValuesData(
  makeRegrTask(data = mldf_wres_br_norm, target = "BrittlenessMPa"), 
  method = c("FSelector_information.gain")
  )

br_fv_plot <- plotFilterValues(regr_fv_br) +
  ggtitle('BR w/ Residual Stress') +
  ylim(0,1) +
  coord_flip()

# No residual
nores_feats = c("Confining Pressure (MPa)", "Peak Stress (MPa)", 
    "Young's Modulus (GPa)", "Poisson's Ratio", 
    "Bulk Density (g/cm3)", "Compressional Sonic (m/s)", "Shear Sonic (m/s)", 
    "Orientation","Group")

mldf_nores_clean <- mldf_nores %>% 
  select(all_of(nores_feats)) %>%
  janitor::clean_names("upper_camel")

mldf_nores_ym_norm = normalizeFeatures(
  mldf_nores_clean %>% select(-PoissonsRatio),
  target ="YoungsModulusGPa")

regr_fv_nores_ym = generateFilterValuesData(
  makeRegrTask(data = mldf_nores_ym_norm, target = "YoungsModulusGPa"), 
  method = c("FSelector_information.gain")
  )

ym_fv_nores_plot <- plotFilterValues(regr_fv_nores_ym) +
  ggtitle('YM w/o Residual Stress') +
  ylim(0,1) + 
  coord_flip()

mldf_nores_pr_norm = normalizeFeatures(
  mldf_nores_clean %>% select(-YoungsModulusGPa),
  target ="PoissonsRatio")

regr_fv_nores_pr = generateFilterValuesData(
  makeRegrTask(data = mldf_nores_pr_norm, target = "PoissonsRatio"), 
  method = c("FSelector_information.gain")
  )

pr_fv_nores_plot <- plotFilterValues(regr_fv_nores_pr) +
  ggtitle('PR w/o Residual Stress') +
  ylim(0,1) +
  coord_flip()

grobs = arrangeGrob(grobs=list(ym_fv_wres_plot,pr_fv_wres_plot,br_fv_plot,ym_fv_nores_plot, pr_fv_nores_plot), widths = unit(c(4,4,4),"in"), heights = unit(c(3,3),"in"))
ggsave("../output/mutual_information.eps", grobs, width = 12, height = 6)
ggsave("../output/mutual_information.pdf", grobs, width = 12, height = 6, dpi=600)

ScottHMcKean/duvernay_geomech_model documentation built on Sept. 12, 2022, 10:08 a.m.