Movatterモバイル変換


[0]ホーム

URL:


Plotting geological/stratigraphicalpatterns

In 2006, theU.S. GeologicalSurvey and theGeologic DataSubcommittee of theFederalGeographic Data Committee established theDigitalCartographic Standard for Geologic Map Symbolization. This is theNational Standard for the digital cartographic representation ofgeologic map features, including line symbols, point symbols, colors,and patterns. The standard patterns that are put forth in this documentare included indeeptime (thanks to thework that Daven Quinn put in to extracting them as SVGs). Let’sexplore the functionality thatdeeptime has to retrieveand utilize these patterns.

Let’s first load the packages we’ll need:

# Load deeptimelibrary(deeptime)# Load rmacrostrat to get datalibrary(rmacrostrat)# Load packages for data visualizationlibrary(ggplot2)library(ggrepel)

Stratigraphic columns

Thermacrostratpackage allows users to access the Macrostrat API, which includesvarious geological data (e.g., lithostratigraphic units) anddefinitions/metadata associated with those data. The package includes anumber ofvignettes thatwalk through how to retrieve and visualize various types of data fromthe database. Here, we’re going to be plotting a stratigraphic columnfor the San Juan Basin, a large structural depression which spans partsof New Mexico, Colorado, Utah, and Arizona. The details aboutdownloading this data are thoroughly presented inthisrmacrostrat vignette. For the purposes of thisdeeptime vignette, we’ll skip ahead and download theunit-level stratigraphic data for this basin during the Cretaceous:

# Using the column ID, retrieve the units in the San Juan Basin columnsan_juan_units<-get_units(column_id =489,interval_name ="Cretaceous")# See the column names and number of rowscolnames(san_juan_units)##  [1] "unit_id"          "section_id"       "col_id"           "project_id"##  [5] "col_area"         "unit_name"        "strat_name_id"    "Mbr"##  [9] "Fm"               "Gp"               "SGp"              "t_age"## [13] "b_age"            "max_thick"        "min_thick"        "outcrop"## [17] "pbdb_collections" "pbdb_occurrences" "lith"             "environ"## [21] "econ"             "measure"          "notes"            "color"## [25] "text_color"       "t_int_id"         "t_int_name"       "t_int_age"## [29] "t_prop"           "units_above"      "b_int_id"         "b_int_name"## [33] "b_int_age"        "b_prop"           "units_below"      "strat_name_long"## [37] "refs"             "clat"             "clng"             "t_plat"## [41] "t_plng"           "b_plat"           "b_plng"nrow(san_juan_units)## [1] 17

Plot stratigraphic column

We now have information for each of the 17 Cretaceous geologic unitscontained within the San Juan Basin, including the age of the top andbottom of each, which is what we will use to plot our stratigraphiccolumn. We’ll follow thermacrostrat vignette to plotthe stratigraphic column:

# Specify x_min and x_max in dataframesan_juan_units$x_min<-0san_juan_units$x_max<-1# Tweak values for overlapping unitssan_juan_units$x_max[10]<-0.5san_juan_units$x_min[11]<-0.5# Add midpoint age for plottingsan_juan_units$m_age<- (san_juan_units$b_age+ san_juan_units$t_age)/2ggplot(san_juan_units,aes(ymin = b_age,ymax = t_age,xmin = x_min,xmax = x_max))+# Plot units, colored by rock typegeom_rect(fill = san_juan_units$color,color ="black")+# Add text labelsgeom_text_repel(aes(x = x_max,y = m_age,label = unit_name),size =3.5,hjust =0,force =2,min.segment.length =0,direction ="y",nudge_x =rep_len(x =c(2,3),length.out =17))+# Add geological time scalecoord_geo(pos ="left",dat =list("stages"),rot =90)+# Reverse direction of y-axisscale_y_reverse(limits =c(145,66),n.breaks =10,name ="Time (Ma)")+# Remove x-axis guide and titlescale_x_continuous(NULL,guide =NULL)+# Choose theme and font sizetheme_classic(base_size =14)+# Make tick labels blacktheme(axis.text.y =element_text(color ="black"))

Use stratigraphic patterns

Isn’t that snazzy! You’ll see that we used the colors that comestraight from Macrostrat (i.e., the “color” column) to visualize thedifferent units. Presumably these colors represent something about thedifferent lithologies of the units. However, we could also representthese lithologies with our standardized FGDC patterns. First, we’ll needto get the lithology definitions from Macrostrat, which includes theFGDC pattern codes:

# Get lithology definitionsliths<-def_lithologies()head(liths)##   lith_id          name          type          group       class   color fill## 1       1 siliciclastic siliciclastic                sedimentary #FFF400  670## 2       2        gravel siliciclastic unconsolidated sedimentary #FFAB00  601## 3       3          sand siliciclastic unconsolidated sedimentary #FFCA00  607## 4       4          silt siliciclastic unconsolidated sedimentary #919AA3  619## 5       5           mud siliciclastic unconsolidated sedimentary #A7ACB0  641## 6       6     claystone siliciclastic       mudrocks sedimentary #B2B5B7  620##   t_units## 1     347## 2    1316## 3   12663## 4    9839## 5    8472## 6    5920

In this output, the “fill” column corresponds to the FGDC patterncode for each lithology. Now, we need to figure out the lithology ofeach of our San Juan units. However, you’ll notice that some units havemultiple lithologies:

# Count lithologies for each unitsapply(san_juan_units$lith, nrow)##  [1] 1 1 1 2 1 1 1 3 3 1 1 1 1 1 1 6 3# Inspect one with multiple rowssan_juan_units$lith[16]## [[1]]##         atts      name   prop lith_id            type       class## 1                shale 0.0714       8   siliciclastic sedimentary## 2            quartzite 0.3571      85 metasedimentary metamorphic## 3            sandstone 0.3571      10   siliciclastic sedimentary## 4        red claystone 0.0714       6   siliciclastic sedimentary## 5 calcareous      clay 0.0714      93   siliciclastic sedimentary## 6 calcareous      clay 0.0714      93   siliciclastic sedimentary

That’s a lot of different lithologies for this single unit! We couldtheoretically represent all of the lithologies, but let’s just go aheadand extract the primary lithology, or the one with the highest “prop”value (i.e., proportion of unit):

# Get the primary lithology for each unitsan_juan_units$lith_prim<-sapply(san_juan_units$lith,function(df) {  df$name[which.max(df$prop)]})table(san_juan_units$lith_prim)#### conglomerate    graywacke  litharenite    quartzite    sandstone        shale##            1            1            1            1            6            7

Now that we’ve assigned a primary lithology to each unit, we canassign an FGDC pattern code to each unit:

# Assign pattern codesan_juan_units$pattern<-factor(liths$fill[match(san_juan_units$lith_prim, liths$name)])table(san_juan_units$pattern,exclude =NULL)#### 602 607 620 654 702##   1   6   7   2   1

Excellent, we now have a pattern code for each unit! Now we can goahead and use these instead of the colors for the fills of the units. Inorder to convert the codes to visual patterns, we’ll need to move thefill argument to a call toaes() and we’llneed to add a call toscale_fill_geopattern():

# Plot with pattern fillsggplot(san_juan_units,aes(ymin = b_age,ymax = t_age,xmin = x_min,xmax = x_max))+# Plot units, patterned by rock typegeom_rect(aes(fill = pattern),color ="black")+scale_fill_geopattern()+# Add text labelsgeom_text_repel(aes(x = x_max,y = m_age,label = unit_name),size =3.5,hjust =0,force =2,min.segment.length =0,direction ="y",nudge_x =rep_len(x =c(2,3),length.out =17))+# Add geological time scalecoord_geo(pos ="left",dat =list("stages"),rot =90)+# Reverse direction of y-axisscale_y_reverse(limits =c(145,66),n.breaks =10,name ="Time (Ma)")+# Remove x-axis guide and titlescale_x_continuous(NULL,guide =NULL)+# Choose theme and font sizetheme_classic(base_size =14)+# Make tick labels blacktheme(axis.text.y =element_text(color ="black"))

Well look at that! Oh, but that legend isn’t very helpful at themoment, is it? Let’s go ahead and customize it to define whatlithologies the different patterns actually represent. We’ll use thecodes and names inliths to set the breaks and labels forthe legend, and we’ll move it to the bottom so it doesn’t mess up ourlabels:

# Plot with pattern fillsggplot(san_juan_units,aes(ymin = b_age,ymax = t_age,xmin = x_min,xmax = x_max))+# Plot units, patterned by rock typegeom_rect(aes(fill = pattern),color ="black")+scale_fill_geopattern(name =NULL,breaks =factor(liths$fill),labels = liths$name)+# Add text labelsgeom_text_repel(aes(x = x_max,y = m_age,label = unit_name),size =3.5,hjust =0,force =2,min.segment.length =0,direction ="y",nudge_x =rep_len(x =c(2,3),length.out =17))+# Add geological time scalecoord_geo(pos ="left",dat =list("stages"),rot =90)+# Reverse direction of y-axisscale_y_reverse(limits =c(145,66),n.breaks =10,name ="Time (Ma)")+# Remove x-axis guide and titlescale_x_continuous(NULL,guide =NULL)+# Choose theme and font sizetheme_classic(base_size =14)+# Make tick labels blacktheme(legend.position ="bottom",legend.key.size =unit(1,'cm'),axis.text.y =element_text(color ="black"))

Note that we also could have usedfgdc_dict() to get alabeling function to use withinscale_fill_geopattern()here, but we want to use the macrostrat lithology names here, not thelithology names from the FGDC standard definitions, so we’ll stick withgrabbing them fromliths.

Further customization

If you usescale_fill_geopattern(), the patterns areused with the default parameters (e.g., line color, background color,and scale). However, if you’d like to customize how the patterns look,you can use thegeom_ functions from theggpatternpackage (e.g.,geom_rect_pattern()) and use the “geo”pattern. Once you have implemented this, any of thefollowing ggplot2 aesthetics can be used to tweak the patterns.Notethat the defaults for these two methods are often different:

aestheticdescriptionscale_fill_geopattern() defaultggpattern default
pattern_typeCode for the FGDC pattern to use“101”“101”
pattern_colorColor used for strokes and pointsoriginal color“grey20”
pattern_fillColor used to fill various closed shapes (e.g., circles) in thepatternoriginal color“grey80”
pattern_alphaAlpha transparency for pattern11
pattern_scaleScale of the pattern (higher values zoom in on pattern)21
fillBackground color“white”“white”

For the pattern code (i.e.,pattern\_type), see the“pattern numbers” in thefullpattern chart for the full list of options. Daven Quinn has alsoassembled more accessible documentation of themappatterns/codes andlithologypatterns/codes.def_lithologies() can also be used tolook up pattern codes for various lithologies (see the “fill” column).Note that codes associated with color variants (e.g., “101-M”) aresupported but will result in the default color variant instead (usuallyblack and white). Most, if not all, color variants can be recreated byadjusting the color, fill, and background of the pattern.

Let’s try increasing the scale of the patterns in our stratigraphiccolumn. We’ll need to switch to usinggeom_rect_pattern().Since we are specifying the patterns within anaes() call,we can usescale_pattern_type_identity() to use those rawpattern codes:

# Plot using ggpatternlibrary(ggpattern)## Warning: package 'ggpattern' was built under R version 4.5.2ggplot(san_juan_units,aes(ymin = b_age,ymax = t_age,xmin = x_min,xmax = x_max))+# Plot units, patterned by rock typegeom_rect_pattern(aes(pattern_type = pattern),pattern ="geo",pattern_scale =4)+scale_pattern_type_identity(name =NULL,guide ="legend",breaks =factor(liths$fill),labels = liths$name)+# Add text labelsgeom_text_repel(aes(x = x_max,y = m_age,label = unit_name),size =3.5,hjust =0,force =2,min.segment.length =0,direction ="y",nudge_x =rep_len(x =c(2,3),length.out =17))+# Add geological time scalecoord_geo(pos ="left",dat =list("stages"),rot =90)+# Reverse direction of y-axisscale_y_reverse(limits =c(145,66),n.breaks =10,name ="Time (Ma)")+# Remove x-axis guide and titlescale_x_continuous(NULL,guide =NULL)+# Choose theme and font sizetheme_classic(base_size =14)+# Make tick labels blacktheme(legend.position ="bottom",legend.key.size =unit(1,'cm'),axis.text.y =element_text(color ="black"))

Oh no! The patterns are more zoomed in as we wanted, but now thecolors are all messed up! This is a good example of how the defaultschange between the two methods. We can easily fix this by specifying ourown defaults:

ggplot(san_juan_units,aes(ymin = b_age,ymax = t_age,xmin = x_min,xmax = x_max))+# Plot units, patterned by rock typegeom_rect_pattern(aes(pattern_type = pattern),pattern ="geo",pattern_color ="black",pattern_fill ="white",fill ="white",pattern_scale =4)+# Use identity of pattern_type aesthetic to set pattern type# Also, substitute lithology names for codes in the legendscale_pattern_type_identity(name =NULL,guide ="legend",breaks =factor(liths$fill),labels = liths$name)+# Add text labelsgeom_text_repel(aes(x = x_max,y = m_age,label = unit_name),size =3.5,hjust =0,force =2,min.segment.length =0,direction ="y",nudge_x =rep_len(x =c(2,3),length.out =17))+# Add geological time scalecoord_geo(pos ="left",dat =list("stages"),rot =90)+# Reverse direction of y-axisscale_y_reverse(limits =c(145,66),n.breaks =10,name ="Time (Ma)")+# Remove x-axis guide and titlescale_x_continuous(NULL,guide =NULL)+# Choose theme and font sizetheme_classic(base_size =14)+# Make tick labels black and increase legend key sizetheme(legend.position ="bottom",legend.key.size =unit(1,'cm'),axis.text.y =element_text(color ="black"))

There we go, now the colors are fixed! Hmm…now it looks a littleboring. Let’s go ahead and bring back those colors from before. We’llmake both the background color (fill) and the pattern fillcolor (pattern_fill) based on the “color” column:

ggplot(san_juan_units,aes(ymin = b_age,ymax = t_age,xmin = x_min,xmax = x_max))+# Plot units, patterned by rock typegeom_rect_pattern(aes(pattern_type = pattern,fill = color,pattern_fill = color),pattern ="geo",pattern_color ="black",pattern_scale =4)+# Need to override the legend defaultsscale_pattern_type_identity(name =NULL,guide =guide_legend(override.aes =list(pattern_fill ="white",fill ="white")                              ),breaks =factor(liths$fill),labels = liths$name)+# Use the raw color valuesscale_fill_identity()+scale_pattern_fill_identity()+# Add text labelsgeom_text_repel(aes(x = x_max,y = m_age,label = unit_name),size =3.5,hjust =0,force =2,min.segment.length =0,direction ="y",nudge_x =rep_len(x =c(2,3),length.out =17))+# Add geological time scalecoord_geo(pos ="left",dat =list("stages"),rot =90)+# Reverse direction of y-axisscale_y_reverse(limits =c(145,66),n.breaks =10,name ="Time (Ma)")+# Remove x-axis guide and titlescale_x_continuous(NULL,guide =NULL)+# Choose theme and font sizetheme_classic(base_size =14)+# Make tick labels black and increase legend key sizetheme(legend.position ="bottom",legend.key.size =unit(1,'cm'),axis.text.y =element_text(color ="black"))

Well isn’t that a nice looking visualization of the stratigraphiccolumn!


[8]ページ先頭

©2009-2025 Movatter.jp