Vostock Ice Core Lab 2020
[Natalie Smith] | [Sophia Vitello]
[Jan 27th 2020]
Analysis of Vostok Ice Core Data
One of the primary results of paleoclimate research over the past decade has been strong evidence for human-influenced (anthropogenic) global climate change. Some of our strongest evidence comes from ice cores taken from undisturbed ice sheets, such as those in Antarctica. Ice cores provide uninterrupted information on important properties of paleoclimate, including local temperature and precipitation rate, humidity, and wind speed as well as changes in atmospheric composition captured at the time that layer of ice froze. In this lab we will work with the original data from the famous Vostok Ice Cores, which carries the distinction of being “the only ice cores that scientists are certain have remained undisturbed for the last interglacial and the penultimate glacial periods” (AGU, 1995).
Vostok is the Antarctic research base founded by the Soviet Union in 1957. Even though the US has a research base at the South Pole, the (now) Russian base at Vostok is much more important scientifically. First of all it sits over one of the thickest (>4km) portions of the East Antarctic Ice Sheet. Second, this ice is relatively undisturbed (elsewhere the ice has been subject to substantial flow and/or seasonal conditions that might destroy the ice). Third, under the 4km ice sheet at Vostok lies a lake the size of Lake Ontario that has been cut off from surface conditions for perhaps 15 million years.
OK - Let’s get started!!
vos<-read.csv("totalvosdataNEW (2).csv")
summary (vos)## depth.m ice_age.ka deltaD.permil dust.cm3g
## Min. : 130.0 Min. : 4.26 Min. :-498.6 Min. : 5.00
## 1st Qu.: 612.5 1st Qu.: 35.57 1st Qu.:-484.1 1st Qu.: 24.91
## Median :1095.0 Median : 73.98 Median :-473.9 Median : 57.00
## Mean :1095.0 Mean : 74.78 Mean :-469.4 Mean :108.21
## 3rd Qu.:1577.5 3rd Qu.:114.05 3rd Qu.:-458.8 3rd Qu.:131.40
## Max. :2060.0 Max. :160.73 Max. :-421.0 Max. :668.43
##
## gas_age.ka CO2.ppmv CH4.ppbv
## Min. : 1.90 Min. :179.9 Min. :318.2
## 1st Qu.: 31.73 1st Qu.:204.7 1st Qu.:422.8
## Median : 70.22 Median :226.3 Median :459.3
## Mean : 71.44 Mean :229.6 Mean :482.8
## 3rd Qu.:110.72 3rd Qu.:254.9 3rd Qu.:546.0
## Max. :156.65 Max. :291.7 Max. :693.4
## NA's :2Load in your data - create a code chunk (+C above - it will appear wherever your cursor is) to load in the vostock ice data and run a summary, the same way you did in the R-markdown Introduction (you can even launch your R-markdown intro html report and copy the code directly into the code chunk you create). Look at the first five rows of the table to make sure it loaded properly using the head() function.
Plotting Information
Okay, now that we’ve loaded in the data, let’s load some libraries that will provide extra functions for us to use to explore the data further! The following code will load the “tidyverse” library which provides MANY useful functions for data science.
As always, when we load a package, we use the function library(). Let’s try those now with the “tidyverse” library.
library(tidyverse) # you don't need quotes when loading the libraryPLEASE TYPE ALL YOUR ANSWERS BETWEEN 2 ASTERISKS SO THEY APPEAR BLUE IN THIS WINDOW AND BOLD IN THE HTML REPORT THAT WILL GET GENERATED AT THE END!!
Question 1: What do you expect to see and why in terms of the relationship between ice age and gas age (as a function of ice depth) based on what you know about sintering (0.5pts)?
Mean ice depth would be a greater value in the ice age than the gas because older snow has been frozen and compacted and sintering has compressed the porous snow
Now lets plot the age of the ice age and the gas age in terms of the ice core depth using some of the functions we just learned about in our data-carpentry tutorial. Run the ggplot library that is part of the tidyverse (it is probably already running since you used it in your visualizing data exercise).
ggplot(data=vos, aes(x=depth.m, y=age.ka, color=key))+
geom_line(aes(y=ice_age.ka, color="Ice Age"))+
geom_line(aes(y=gas_age.ka, color="Gas Age"))+
ggtitle("Age of Ice and Gas by Ice Depth")
Now we should have a nice figure of ice and gas ages by ice core sampling depth. Let’s answer a few questions:
Question Series 2.
2a) How much younger, in years, is a bubble of gas estimated to be compared to the ice that surrounds it at a depth of 1000 meters? Hint: you can insert a text chunk here, then use the code “filter(vos, depth.m==1000).” The filter function allows us to see the data at specific values (0.5 pt).
filter(vos, depth.m==1000)## depth.m ice_age.ka deltaD.permil dust.cm3g gas_age.ka CO2.ppmv CH4.ppbv
## 1 1000 66.01 -483.4 25 61.86 192.84 425.414150 years
2b) How much variability is in the difference between ice and gas ages over depth? Hint: you can use the code “mutate(vos, diff_age=ice_age.ka-gas_age.ka).” The mutate function allows us to add new columns to our data where we can then calculate values from existing columns. Use the “summary()” function to look at the variabilty in this new “diff_age” variable over depth (0.5 pts).
vos<-mutate(vos, diff_age=ice_age.ka-gas_age.ka)
summary(vos$diff_age)## Min. 1st Qu. Median Mean 3rd Qu. Max.
## 2.230 2.848 3.435 3.336 3.790 4.250Between 4,250 and 2,230 years
2c) Why might age gaps be so large and differ so much (1 pts)?
Rate of sintering, at which snow compacts into ice, is different based on the physical state.
The Temperature Record
In the text block below, you will see html tags to indicate that the isotope values should be super-scripts. This is kind of clunky-looking inside R, but will present very nicely in your final document! Being able to read the html tags will be a useful skill since html is part of so many platforms. Also, the tags produce nicely formatted equations that will actually show if you hover your cursor over them.
Oxygen has two stable isotopes of importance, 16O and 18O, which vary in mass. Differences in the amounts of these isotopes in an oxygen-bearing molecule are measured using a ratio of 18O/16O which is then compared to the ratio of 18O/16O in average ocean water. This comparison is called \(\delta\) 18O (pronounced delta-18-O). Variations in the \(\delta\) 18O of the oxygen in the water molecule, H2O, can be useful in understanding the hydrological cycle. Average ocean water has a value of 0‰ \(\delta\) 18O (‰ is pronounced permil and is the symbol for one-thousandth. It is analogous to %, percent, which is the symbol for one-hundredth). See Figure 1: https://georgetown.instructure.com/courses/102116/pages/vostock-lab-figures (when working inside R, you’ll need to copy and paste the link into your browser, or you can [shift-mouseclick], once your html report is generated, this will appear as a functioning hyper-link).
When ocean water evaporates, water with the lighter oxygen isotope (16O) evaporates more easily because it is lighter than a water molecule with the heavier oxygen isotope (18O). Therefore, water vapor in the atmosphere ends up with a smaller percentage of 18O than ocean water. Its \(\delta\) 18O is a few permil negative, say around –3‰. When that vapor passes over land and condenses to form rain, the heavier isotopes that did make it into the clouds condense to a greater extent than the lighter isotopes, again due to mass. Precipitation that falls then has a more negative \(\delta\) 18O than seawater, but more positive \(\delta\) 18 than the clouds from which it fell. This effect is more pronounced in cold climates than in warmer ones, because temperature is what drives the evaporation/condensation processes. Snow has much more negative \(\delta\) 18O than rain. Since temperature drives this process, an equation has been derived to relate temperature to \(\delta\) 18O.
The scale we use here to describe the variation in isotopic ratios of substances is \(\delta\) D (pronounced delta-D). For most precipitation on Earth, the following relationship applies:
\(\delta D = 8 * \delta ^{18}O + 9\)
We want to convert the isotopic ratios to temperature changes (delta temperatures – this means “change in” temperature, it is not the absolute temperature) that describe variations in the temperature of the ocean from which the ice was originally evaporated. The formula to do this is:
\(\delta Temperature = (\delta D + 440)/6.2\)
Let’s add another column to the “vos” data for the calculation of Delta.T, you can do it the same way we created the “diff_age” column a few minutes ago by using the “mutate()” function, only now call it Delta.T. As always, its a good idea to check and make sure the table generated with your new field. There are lots of ways you can do that (click on table name in above right pane, do summary() function or head() function are a few)
vos<-mutate(vos,Delta.T=(deltaD.permil+440)/6.2)
head(vos)## depth.m ice_age.ka deltaD.permil dust.cm3g gas_age.ka CO2.ppmv CH4.ppbv
## 1 130 4.26 -438.94 6.58 1.90 274.00 NA
## 2 140 4.66 -435.94 12.42 2.28 273.07 NA
## 3 150 5.06 -439.72 11.00 2.66 272.14 666.93
## 4 160 5.45 -441.88 7.47 3.04 271.22 653.54
## 5 170 5.85 -434.77 6.13 3.41 270.29 640.15
## 6 180 6.27 -440.36 11.57 3.82 268.39 594.65
## diff_age Delta.T
## 1 2.36 0.17096774
## 2 2.38 0.65483871
## 3 2.40 0.04516129
## 4 2.41 -0.30322581
## 5 2.44 0.84354839
## 6 2.45 -0.05806452vos<-mutate(vos,Temp.degC=(Delta.T-55.5))
summary(vos$Temp.degC)## Min. 1st Qu. Median Mean 3rd Qu. Max.
## -64.96 -62.60 -60.98 -60.24 -58.54 -52.43Now that we’ve created that new variable Delta.T, let’s plot the temperature curve as a function of ice depth (or age).
plot(vos$ice_age.ka, vos$Delta.T, type ="l", main = "Temperature Change by Age", xlab = "Ice Age (K)", ylab = "Estimated change in temperature C", xaxt="n")
axis(1, at = seq(0, 150, by = 10), las = 2)
Question Series 3
3a) How many degrees has the climate varied in the past 150,000 years (0.5 pt)?
summary(vos$"Delta.T")## Min. 1st Qu. Median Mean 3rd Qu. Max.
## -9.456 -7.105 -5.476 -4.736 -3.039 3.06812.524
3b) According to the graph, approximately when was the last period of maximum ice cover (glacial maximum) (0.5 pt)?
20,000 years ago
3c) When did the last interglacial time period begin (the time that the temperature was on average warmer than today) (0.5pts)?
10,000 years ago
What Can We Learn about Climate Today?
Let’s add another column to our dataset (using the mutate() function) to estimate the temperature at Vostok. To do this subtract 55.5 degrees from the numbers in the delta temperature column (Delta.T) to get an estimate of the Vostok temperature itself. Call this variable “Temp.degC” and check your table to make sure it added properly
Now we’ll create scatterplots of CO2 (x-axis) vs. temp.degC (y-axis) and CH4 (x-axis) vs. temp.degC (y-axis) using the Vostok data. We will also be adding in trendlines using the abline() and lm() functions. The lm() function determines what the best fit line is (and will also tell us if it is “statistically significant” and also how much variability it explains) and the abline() function just plots the best fit line. Don’t worry too much about these details right now!
par(mfrow=c(1,2)) # this lets us show multiple plots at once
plot(vos$CO2.ppmv, vos$Temp.degC) # plots the points
abline(lm(vos$Temp.degC~vos$CO2.ppmv)) # plots the best fit line
plot(vos$CH4.ppbv, vos$Temp.degC) # plots the points
abline(lm(vos$Temp.degC~vos$CH4.ppbv)) # plots the best fit line
Next, we will analyze the goodness-of-fit (R2) of our bivariate plots using linear regression. The R-squared value (R2) is a statistical measure of how well a regression line approximates real data points; an R-squared of 1.0 (100%) indicates a perfect fit. So if you had an R2 of 1.0 for the CO2 vs. temperature plot, 100% of the variation in temperature could be explained by knowing the CO2 concentration.
cotwo <- lm(vos$Temp.degC~vos$CO2.ppmv) # reruns the linear regression and puts the results into an object called “cotwo”
summary(cotwo) # summarizes the results of your test##
## Call:
## lm(formula = vos$Temp.degC ~ vos$CO2.ppmv)
##
## Residuals:
## Min 1Q Median 3Q Max
## -5.6186 -1.1265 0.1969 1.2913 4.2469
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -81.045093 1.093529 -74.11 <2e-16 ***
## vos$CO2.ppmv 0.090634 0.004728 19.17 <2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 1.83 on 192 degrees of freedom
## Multiple R-squared: 0.6568, Adjusted R-squared: 0.655
## F-statistic: 367.4 on 1 and 192 DF, p-value: < 2.2e-16chfour <- lm(vos$Temp.degC~vos$CH4.ppbv) # runs the linear regression and puts the results into an object called “cofour”
summary(chfour) # summarizes the results of your test##
## Call:
## lm(formula = vos$Temp.degC ~ vos$CH4.ppbv)
##
## Residuals:
## Min 1Q Median 3Q Max
## -7.5753 -1.1699 -0.1495 1.4518 4.3674
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -73.680919 0.893139 -82.50 <2e-16 ***
## vos$CH4.ppbv 0.027737 0.001824 15.21 <2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 2.078 on 190 degrees of freedom
## (2 observations deleted due to missingness)
## Multiple R-squared: 0.5491, Adjusted R-squared: 0.5467
## F-statistic: 231.3 on 1 and 190 DF, p-value: < 2.2e-16Question Series 4
4a) What is the Adjusted R2 value for CO2 (0.25 pts)?
.655
4b) What is the Adjusted R2 value for CH4 (0.25 pts)?
.5467
4c) Compare these two values. How are they different and what does that tell you about the relationship between greenhouse gases and temperature (1 pt)?
The R squared value for Temperature vs. CO2 in the atmosphere is greater, suggesting CO2 to be the primary greenhouse gas, compared to methane, influencing changes in temperature of the atmosphere
4d) Why are the values so much lower than 1.0 (1 pt)?
The values are so much lower than 1.0 because there is a wide breadth of factors contributing to atmospheric temperature, including additional greenhouse gases
Comparing to Today’s CO2 Concentration
Question Series 5
5a) Predict the temperature at Vostok today. First, look up and report today’s CO2 concentration (0.5 pts):
Today’s CO2 concentration is: 414.67ppmv
Hint: If today’s value is not available, you can scroll down and it should show you the value from a couple of days ago.
Use the results from your linear regression to estimate the current temperature at Vostock from the relationship between CO2 and temperature that you modeled above. Your model is the line you drew in the scatter plot. Remember that y = mx + b where “m” is the slope and “b” is the intercept.
5b) What do y, x m, and b represent in terms of this specific example (0.5 pts ea)?
Y is Temperature that we calculate -43.522 X is CO2 concentration that I looked up 414.67ppmv m is is the average rate of temperature change as temperature increases .091 b is at 0 CO2 concentration the value of the temperature -81.05
Remember, you have already output the results of your lm() function - and that will also tell you the estimates for the parameters of your model (in other words, the value of “m”, the slope and the value of “b”, the intercept).
5c) Based on the equation y = mx +b, what do you predict is today’s temperature at the Vostock Ice station in Antarctica (hint, you can create a chunk of code and just calculate it as an arithmetic function) (1 pt)?
(.0906*414)-81.045## [1] -43.53665d) What is today’s ACTUAL temperature at Vostock in degrees centigrade (0.5 pts). You can look this up at: https://www.accuweather.com/en/aq/vostok-station/2273742/weather-forecast/2273742)?
-31.11
5e) Give two reasons why the observed value might be different than your prediction (0.5 pt ea)?
The temperature can vary based on the time of day, different weather conditions such as wind, precipitation, cloud patterns
5f) You can get a more complete view of Vostock’s current climate here: https://en.wikipedia.org/wiki/Vostok_Station#Climate. Of the available temperatures that you could use to compare your prediction, which would you choose and why do you think is the most fair comparison (1 pt)?
I chose the daily average for the month of February because it is the end of the month and the daily average would represent the most accurate temperature calculation of this time period
5g) What is the temperature value you chose for comparison and how did it compare to your prediction? (0.5 pts)
-44.3 degrees celsius - very close value with a difference of .8 degress celsius
OK - YOU’RE DONE WITH YOUR FIRST R-MARKDOWN EXERCISE. IF EVERYTHING WORKED WHEN YOU RAN YOUR CODE ABOVE, YOU SHOULD BE ABLE TO “KNIT” THIS SCRIPT AND ASSOCIATED OUTPUT INTO AN HTML DOCUMENT. FIRST, SAVE YOUR DOCUMENT (JUST HIT THE LITTLE DISK ICON) THEN CLICK ON “KNIT” ON THE TASK BAR ABOVE (3 pts for correctly inserting and implementing code and outputting this report)