#  Dynamic Data in the Statistics Classroom

Johanna Hardin

##  1  Introduction

There has been a recent push to change the way we - as statisticians - engage
pedagogically with complex real-world data analysis problems. Two parallel
forces have directed us toward embracing more complex-data real-world problems
in the classroom.

The first driver is a call from educators to make statistics more relevant to
students’ experiences. In order to address the perspective of students who
previously have been reported to “exhibit remarkably little curiosity about
the material they are analyzing", Brown and Kass (2009) suggest the importance
of “real-world problem solving" to get students engaged in their analyses.
Gould (2010) wants introductory statistics students to leave the course with
“a set of … attitudes about data that are immediately applicable to their
lives." Gould and Çetinkaya-Rundel (2013) suggest putting “data at the center
of the curriculum.” Horton _et al._ (2015) assert that “statistics students
need to develop the capacity to make sense of the staggering amount of
information collected in our increasingly data-centered world.” And Zhu _et
al._ (2013) remind us that “data pre-processing bridges the gap from data
acquisition to statistical analysis but has not been championed as a relevant
component in statistics curricula.”

The second driver arises from students and other stakeholders and is harder to
document with references. Certainly my own experience (and that of my
colleagues) is that students engage at the deepest level when they (a) care
about the problem, and (b) understand the data collection, study design, or
motivation of the problem. Kuiper and Sturdivant (2015) report, “Our student
evaluations of these materials support Gould’s (2010) comments suggesting that
previously collected and cleaned data were considered abstract to the student
… We have found that unless students 1) have collected the data themselves or
2) clearly see where and how the data were collected, they often fail to
appreciate it."

Indeed, even if the work is done for them, as Grimshaw (2015) argues, there is
much benefit in the students seeing how the data were procured:

> Using the vocabulary of Wickham (2014), teachers hide the ‘messy data’
> aspects and provide ‘tidy data’ - even when students possess the data skills
> required to work with the messy data. It is valuable for students to not
> only have many authentic data experiences but also to have the professor
> model the correct application of statistics by showing work with messy data
> in lectures…What is good for statistics majors can also be applied to
> introductory courses. There may be no data skills on the learning outcomes
> for these courses, but some examples and homework in an introductory course
> may be modified and/or updated to use the original source data instead of a
> curated dataset. The data skills required would certainly be modest and need
> to fit student backgrounds. The objective would be for students to see that
> data skills are required in an analysis. Students may rely on code provided
> to them that results in their own copy of the dataset.

There have been repeated calls to use real data in the classroom (Cobb, 1991,
1992, 2007, 2011; Workgroup on Undergraduate Statistics, 2000; American
Statistical Association Undergraduate Guidelines Workgroup, 2014; GAISE
College Group, 2005; Carver _et al._ , 2016). There are benefits to using real
data sets (as opposed to made up numbers) in textbooks and in the classroom.
Indeed, there are virtually no textbooks of any kind (AP Statistics,
Introductory Statistics, second course in statistics) being written today
without the vast majority of examples taken from actual studies or databases.
Additionally, there has been a push to infuse R with datasets that are
relevant, recent, and sophisticated (e.g., nycflights13 and other mosaic
datasets, (Pruim _et al._ , 2014)).

Unfortunately, by nature, however, all of the data given in a textbook or an R
package is static. Currently, there do not exist mechanisms to continually
update any dataset provided by the course materials. The most comprehensive
baseball dataset compiled and provided to the students will be out of date
(and less interesting) by next fall. The good news, however, is that, outside
of the statistics classroom, data of all kinds are being updated in real time
publicly and accessibly. And even better news is that R developers are
continually improving interfaces to the vast amounts of public data. For
example, the `tidyverse` package imports a handful of packages that make
downloading data straightforward (e.g., see Wickham (2014) and more recently
the packages: `rvest`, `readr`, `readxl`, `haven`, `httr`, and `xml2`; expect
more to come).

Grimshaw (2015) delineates data along two axes: the first axis describes the
source and format for the data; the second axis describes the amount of
wrangling required on the data. Both axes are scored as good/better/best.
Indeed, Grimshaw gives examples, sources, and R code for several examples.
However, he included only one dataset from the “best/best" category that
worked well in his classroom. Out of the three “best/best" datasets he used,
two were not well received because of excessive time, computing expertise, or
domain knowledge needed to wrangle the data into usable format and extract
meaning from the dataset.

In this manuscript, we seek to address the difficulty of using dynamic data in
the classroom by curating additional best/best datasets (including full R
Markdown files to facilitate reproducible analysis). The available resources
aim to open the world of dynamic data to those who have not previously worked
directly with downloaded data. We hope our work can act as a starting point
for those interested in building their data scraping skills. Additionally, the
R Markdown files can be used as scaffolding for assignments designed to help
students engage with online resources.

In the next section we work through an entire example including showing the
steps and R code for downloading the data, providing example code for using it
in a typical introductory statistics classroom, and suggesting ideas for
expanded use beyond introductory statistics. Section 3 gives summaries of
additional dynamic datasets that are fully curated as R Markdown files in the
supplementary resources on Github (https://github.com/hardin47/DynamicData) as
well as a list of other places to look for dynamic data, and we provide
closing thoughts and ideas for future work in the Conclusion. The Appendix
describes a few technical details that are helpful for making dynamic data
accessible to students at all levels.

##  2 Complete Example: College Scorecard Data

Supplementary materials provide full R Markdown files for complete analysis on
five different dynamic datasets. We start with a complete description of the
College Scorecard materials as a way of illustrating the available resources.

Data on characteristics of US institutions of higher education were collected
in an effort to make more transparent issues of cost, debt, completion rates,
and post-graduation earning potential. An undertaking of the U.S. Department
of Education, the College Scorecard data represent a compilation of
institutional reporting, federal financial aid reports, and tax information.
The process of gathering and compiling the data is well documented on the
College Scorecard website https://collegescorecard.ed.gov/data/documentation/.
One caveat is that some of the variables have only been collected on students
receiving federal financial aid. Biases inherent to analyses done on data
collected from a subgroup should be considered.

The College Scorecard dataset is incredibly rich. The individual institutions
are broken down by total share of enrollment of various races, family income
levels, first generation status, age of students, etc. Additionally, the data
give matriculation information like SAT scores as well as graduation
information like completion rate and income level. The dataset allows for a
student to investigate political or personal hypotheses about college
education and the costs and benefits within. The variables are described in a
data dictionary given at
https://collegescorecard.ed.gov/assets/CollegeScorecardDataDictionary.xlsx.

###  2.1 Downloading the Data

For each of the five fully curated dynamic datasets, there is an R Markdown
file (available at https://github.com/hardin47/DynamicData) which scrapes the
data from an outside web source (presumably kept current and public by some
other organization). The downloading step is shown here in the manuscript only
for the College Scorecard data. Note, as discussed above, the data are
downloaded directly from the website managed by the US Department of Education
(and not stored locally as a csv file). It is worth pointing out that the code
given here is more complicated than what is standard in an introductory
statistics course (especially at the beginning of the semester!). However,
both the GapMinder and Wikipedia examples can be adjusted in very simple ways,
allowing each student to work with their own dataset. That is, you can
scaffold an assignment by providing the majority of the downloading code and
having the student fill in the URL (Wikipedia example) or variable names
(GapMinder example). The R code for downloading the data for all of the
examples are given in the supplementary R Markdown files, and each uses
slightly different functions and syntax.

First, load the data into R.

    college_url <-
      "https://s3.amazonaws.com/ed-college-choice-public/Most+Recent+Cohorts+(All+Data+Elements).csv"

    college_data <- readr::read_csv(college_url)
    dim(college_data)

    ## [1] 7804 1728

Next, use data wrangling methods to clean and organize the variables.

    college_debt = college_data %>%
      dplyr::select(region, HBCU, DEBT_MDN, md_earn_wne_p10) %>%
      mutate(DEBT_MDN = readr::parse_number(DEBT_MDN),
             md_earn_wne_p10 = readr::parse_number(md_earn_wne_p10)) %>%
      mutate(HBCU = ifelse(HBCU=="NULL", NA, HBCU)) %>%
      mutate(region2 = ifelse(region=="0", "Military",
                      ifelse(region=="1", "New England",
                      ifelse(region=="2", "Mid East",
                      ifelse(region=="3", "Great Lakes",
                      ifelse(region=="4", "Plains",
                      ifelse(region=="5", "Southeast",
                      ifelse(region=="6", "Southwest",
                      ifelse(region=="7", "Rocky Mnts",
                      ifelse(region=="8", "Far West", "Outlying"))))))))))

    summary(college_debt)

    ##      region          HBCU              DEBT_MDN      md_earn_wne_p10
    ##  Min.   :0.000   Length:7804        Min.   :   333   Min.   :  8400
    ##  1st Qu.:3.000   Class :character   1st Qu.:  7710   1st Qu.: 24200
    ##  Median :5.000   Mode  :character   Median :  9833   Median : 31200
    ##  Mean   :4.621                      Mean   : 11830   Mean   : 33233
    ##  3rd Qu.:6.000                      3rd Qu.: 15462   3rd Qu.: 39200
    ##  Max.   :9.000                      Max.   :131335   Max.   :250000
    ##                                     NA’s   :1163     NA’s   :2168
    ##    region2
    ##  Length:7804
    ##  Class :character
    ##  Mode  :character
    ##
    ##
    ##
    ##

###  2.2 Using dynamic data within a typical classroom

Using the downloaded data, we start by applying a technique from the
introductory curriculum to a research question of interest based on the
College Scorecard data. College debt is of particular interest to many college
students, but debt can be mediated by post-graduation income. To fully
investigate the relationship between the variables, we provide both confidence
and prediction intervals for both variables.

After calculating a few individual intervals, we show all intervals
represented graphically and broken down by geographic region. Note that the
visual representations do not represent a simple summary plot of the data, and
we leave it open to the instructor to have the students engage more deeply
with the many available variables.

Using the two variables measuring amount of debt of a typical (i.e., median)
college graduate and median earning 10 years after matriculation, we create
both confidence intervals and prediction intervals – keeping in mind that the
observational unit is an academic institution. Note that the calculations
below are for both confidence and prediction intervals. The confidence
interval agglomerates institutions over the entire dataset; however, the
prediction value is for a single institution (which is the observational
unit). The analysis lends itself nicely to a conversation about confidence vs.
prediction intervals as well as observational units as institution vs. as
individual student. It is worth pointing out to the students that the
prediction intervals likely hold more information related to their individual
experiences than the confidence intervals. However, the unit of prediction is
for an institution, and so the individual student debt and income is likely
even more variable than shown here. Additionally, Figure 1 demonstrates the
effect of samples size: consider the comparison of the Military intervals (one
school) to the intervals for all of the US institutions (about 6000 schools).
(Note: the intervals given in Figure 1 were created using an ANOVA model where
the within variance is calculated across all regions, which is how the
interval for military schools can be calculated. You may or may not want to
bring that up with your students.)

The following R code uses the `mosaic` package to directly calculate both
prediction and confidence intervals. Note the formula interface given by the
tilde is described in detail here:
http://rpruim.github.io/eCOTS2014/Workshop/Modeling.html.

    debt_mod <- lm(DEBT_MDN~1, data = college_debt)
    debt_fun <- mosaic::makeFun(debt_mod)
    debt_fun()

    ##        1
    ## 11829.78

    debt_fun(interval="confidence")

    ##        fit      lwr      upr
    ## 1 11829.78 11692.44 11967.13

    debt_fun(interval="prediction")

    ##        fit      lwr      upr
    ## 1 11829.78 636.2383 23023.33

    earn_mod <- lm(md_earn_wne_p10~1, data = college_debt)
    earn_fun <- mosaic::makeFun(earn_mod)
    earn_fun()

    ##        1
    ## 33232.59

    earn_fun(interval="confidence")

    ##        fit      lwr      upr
    ## 1 33232.59 32864.78 33600.41

    earn_fun(interval="prediction")

    ##        fit      lwr      upr
    ## 1 33232.59 5616.893 60848.29

The intervals are interesting, but they might be even more interesting if
broken down by region and shown visually. Note how much smaller the confidence
intervals are from the prediction intervals! The difference indicates lots of
variability across institutions and large sample sizes.

Figure 1:  The x-axis represents the region of the institution. The y-axis
represents either the amount of debt 10 years after matriculation (orange) or
the amount of income 10 years after matriculation (blue). Confidence intervals
for the average values (within region) are given by the solid lines.
Prediction intervals for individual institutions are given by the dashed
lines. The solid dot represents the center of both types of intervals (broken
down by debt and income).

###  2.3 Thinking outside the box

Additionally, for each of the supplementary dynamic data R Markdown files, we
add a section for each analysis which is based on topics that are not
traditionally taught in introductory statistics classes. The additional
analysis is done not only to expand the tool box of the students but also to
teach the students that they can often think about the problem in
sophisticated ways even if all their tools come only from the introductory
course.

The College Scorecard dataset is incredibly rich and can be used for many
different types of model building: linear, logistic, machine learning. Indeed,
thinking about interaction terms could be particularly insightful. Here, we
give an example of regressing earnings on debt with the interaction term as
whether or not the institution is one of the Historically Black Colleges and
Universities (HBCU). Figure 2 displays the separate regression lines for the
two distinct types of institutions.

Figure 2:  Median income regressed on debt. For the analysis, HBCU is
interacted with debt to provide two distinct (and not parallel) regression
lines. HBCU institutions are given in blue, and non-HBCU institutions are
given in orange.

Many interesting conversations can ensue based on the regression of income on
debt. Reminding the students that each observation is an institution is an
important starting point. Additionally, students should be able to volunteer
the dangers of using a model like this to suggest causality. Last, there might
be room to discuss an inferential analysis of whether HBCUs are statistically
different from non-HBCUs (noting the substantial differences in sample sizes).

It is not hard to come up with additional questions to investigate with the
College Scorecard data. Indeed, because the data relate directly to college
students, they should be able to find ways to engage with the data. We
recommend continued conversations about how the data are valuable to the
larger community, but that the information is not always complete (e.g., many
variables are collected only on students who fill out financial aid forms) and
not causative.

##  3 Dynamic Data Projects

For each of four additional dynamic datasets, we describe the source of the
data and relevant variables & research questions, some standard and graphical
techniques (“using dynamic data within a typical introductory statistics
classroom"), and a statistical analysis appropriate for a course after
introductory statistics (“thinking outside the box"). We also provide source
information for an additional nine dynamic datasets.

###  3.1 Wikipedia Data

Wikipedia stores most of its tabular data in HTML tables. To scrape HTML
tables from any website (or HTML file), use the R function XML::readHTMLTable.
As an initial foray into downloading data directly from the internet into R,
Wikipedia tables provide a nice introduction. In the supplementary R Markdown
file associated with the Wikipedia data analysis, we also walk through some of
the useful aspects of using `dplyr` to wrangle the data111Original idea for
this example provided by Nick Horton, Amherst College..

Figure 3:  The log of the retail value of music sales broken down by whether
or not the sales have increased or decreased (presumably over the previous
year, although the Wikipedia documentation does not specify the time period
over which the change is measured).

#### Using dynamic data within a typical classroom

The Wikipedia analysis given in the fully curated files explores an HTML table
on sales of music (physical and digital) in 2014,
https://en.wikipedia.org/wiki/Music_industry. One variable gives an indication
of whether the retail value of the music sales has increased or decreased.
Using the country-level music data, we perform a t-test, a Wilcoxon rank sum
test, data transformations, and boxplots to investigate music retail sales
(analysis given in supplementary materials, not shown here). By grouping the
data into two categories we can investigate whether there is any statistical
difference between the total average retail sales (in US$) between those
countries for whom retail sales increased versus those that decreased. The
p-value for the initial t-test is reasonably large, but the boxplot shows that
the difference in variability across the two groups is also large with a
sample that either has large outliers or a long skewed right tail. Because the
technical assumptions do not appear to be met, a log transformation of the
data or a non-parametric test might be better assessments of the data (see
Figure 3). The analysis leads to conversations about the source of the data
and the reasons why p-values are non-significant. The example extends easily
to each student choosing their own Wikipedia page and data table, graphical
representations, and statistical analyses.

#### Thinking outside the box

Among the variables in the Wikipedia music dataset are the breakdown
(percentages) of how the retail sales are distributed across physical,
digital, performance rights, and synchronization. We might want to see whether
there is a dependency of total retail sales on the breakdown of types of
products. The problem is not well suited to introductory statistics as there
is not an obvious statistic we can use within a sampling distribution (to
create a p-value, etc.). Because there does not seem to be an obvious
mechanism for evaluating the breakdown of products (and how “different" they
are), we consider an ad-hoc measure and perform a permutation test to assess
significance. The average breakdown of retail sales is given in the R Markdown
file in the supplementary materials. One way to measure a discrepancy between
the retail sales and the consistency of product breakdown is to correlate the
retail sales with the sum of squared distances from the average breakdown of
product types. We see that the metric we created to find a relationship
between retail sales and breakdown of types of product does not show
significance. A student project could be to think about different ways to
measure how the breakdowns can be considered to be different.

###  3.2 NHANES Data

The NHANES data come from the National Health and Nutrition Examination
Survey, surveys given nationwide by the Center for Disease Controls (CDC). The
CDC adopted the following sampling procedure:

  1. 1.

Selection of primary sampling units (PSUs), which are counties or small groups
of contiguous counties.

  2. 2.

Selection of segments within PSUs that constitute a block or group of blocks
containing a cluster of households.

  3. 3.

Selection of specific households within segments.

  4. 4.

Selection of individuals within a household

About 12,000 persons per 2-year cycle were asked to participate in NHANES.
Response rates varied by year, but an average of 10,500 persons out of the
initial 12,000 agreed to complete a household interview. Of these, about
10,000 then participated in data collection at the mobile exam center. The
persons (observational units) are located in counties across the country.
About 30 selected counties were visited during a 2-year survey cycle out of
approximately 3,000 counties in the United States. Each of the four regions of
the United States and metropolitan and non-metropolitan areas are represented
each year. As such, the data collection is ongoing, and the data are updated
on the NHANES website periodically, http://www.cdc.gov/nchs/nhanes.htm. This
manuscript uses the 2011-2012 NHANES data, but we expect the data to be
updated regularly, and the URL should simply change to 2013-2014 when it
becomes available.

Note that the NHANES data are available in the `mosaic` package (Pruim _et
al._ , 2014), but the `mosaic` version of the NHANES data is static (from
2011-2012), and the data has been cleaned with pre-selected variables.
Additionally, the variables can be downloaded directly using the `nhanesA`
package in R.
https://cran.r-project.org/web/packages/nhanesA/vignettes/Introducing_nhanesA.html.
By accessing the data directly from the CDC’s website, students become more
involved in the data analysis process, understanding what they can and cannot
get from the data.

The variables in the CDC’s online NHANES dataset are virtually limitless. We
use a few different datasets, merging them based on an individual identifier
in the dataset. The variable information is all given online, but each at a
different webpage. For example, the demographic data is at
http://wwwn.cdc.gov/nchs/nhanes/search/variablelist.aspx?Component=Demographics&CycleBeginYear=2011.

A further important aspect to the example is that (as described above) the
data do not constitute a simple random sample (it is a weighted sample) from a
population. The sampling scheme can be part of the statistical inquiry into
the data analysis, or it can be set by the instructor in the template used to
download the data. The data which is directly downloaded from the CDC website
includes variables on the weighting scheme. In the R Markdown file provided
with this manuscript, we demonstrate how to create a dataset which can act as
a simple random sample from the population.222The weighting analysis was
motivated by work done by Shonda Kuiper
(http://web.grinnell.edu/individuals/kuipers/stat2labs/weights.html) as well
as the Project Mosaic Team
(https://cran.r-project.org/web/packages/NHANES/NHANES.pdf). The figure and
analysis below are done with a proxy simple random sample.

#### Using dynamic data within a typical classroom

In the supplementary materials (not shown here), we start with a comparison of
body mass index (BMI) for those in committed relationships and those not in
committed relationships. The graphs of the two BMI distributions look quite
similar, and the t-test shows a non-significant difference in means. The
results prompt a discussion about averages versus individual results,
causation, and sample size.

#### Thinking outside the box

Because of the large sample size and the ability to determine the functional
form of a non-linear relationship, smoothing techniques can be used to model
quantitative variables. Adding a smooth curve to a standard scatterplot leads
(see Figure 4) to discussions about how smooth curves are estimated, SE of the
smooth curve, extra variability and instability due to extremes and fewer data
points on the ends. However, extrapolation (note that the two curves have
different ranges?) and the slopes of the two curves not seeming to be
different (no interaction?) might warrant further study.

Figure 4:  A scatterplot of weight and height broken down by gender. A smooth
regression is fit to the points, and the confidence interval for the smooth
curve as well as the increase in variability for small and large values of
height are part of an important discussion in a second level regression
course.

###  3.3 Global Literacy Rate Data

The next dynamic dataset comes from GapMinder http://www.gapminder.org/.
GapMinder also has a plethora of variables from which students can choose
(according to their own interests), but here we work with literacy rates
measured at the country level. The analysis (and more importantly, the data
scraping from GapMinder) could easily be extended for students interested in
all sorts of political, social, environmental, or demographic data available.
Understanding political and demographic trends across both time and location
can provide very interesting insight into economic or political science
questions. Alternatively the GapMinder data can be perused in a descriptive or
graphical manner.

#### Using dynamic data within a typical classroom

The introductory analysis considers gender differences in literacy rates and
uses a linear model on the difference between female and male literacy rates
across time (the analysis is available in the supplementary materials and not
shown here). We show a graphical representation and discuss model assumptions
including sampling and independence of residuals. The model indicates that the
difference between male and female literacy rates is shrinking over time.
However, we worry about the effects of other variables and encourage a more
complete analysis. Indeed, there may be large biases in our model if important
explanatory variables have been left out.

The data provided are ideal for a fantastic classroom conversation about
causation, causal mechanisms, and confounding variables.

#### Thinking outside the box

In the second analysis, we work with the additional variable continent. The
trends observed in the first analysis hold up in the second analysis (i.e.,
the difference declines over time in each of the continents). However, there
are additional considerations to be made, for example, the differences between
the slopes across the continents (the analysis is available in the
supplementary materials and not shown here). We suggest additional
explorations into the independence of the residuals and more advanced spatio-
temporal patterns of literacy.

###  3.4 Weather Data from NOAA Buoys

The National Oceanic and Atmospheric Administration (NOAA) is the American
federal agency in charge of collecting information and making decisions
related to the oceans and the atmosphere. Throughout North America, they
supply weather stations which are located both along the coast as well as in
the middle of the ocean (on buoys). Among other variables, the weather
stations collect information on wind, humidity, temperature, visibility, and
atmospheric pressure. The data is all publicly available on NOAA’s website,
http://www.ndbc.noaa.gov/.

#### Using dynamic data within a typical classroom

Although the data do not constitute a random sample, they are very likely to
be quite representative with respect to the difference in wind and air
temperature at that location over the year. In the supplementary files (not
shown here), we used a paired analysis (i.e., subtract the two variables and
treat them as a single variable) to find a confidence interval for the true
difference in temperature between wind and air. Also, we find a prediction
interval for the difference in temperatures across individual measurements.

#### Thinking outside the box

Although a full analysis of the data would warrant multiple years of data (so
as to understand yearly trends), we can estimate the spectral density of the
time series using a smoothed periodogram (the data below represent
measurements every hour for all of 2014). In the smoothed periodogram (see
Figure 5) the x-axis is the frequency (one over the period) and y-axis
represents the correlation (normalized) between the cosine wave at that
frequency and the time series. We can see that wind speed has strong
correlation at period 12 hours and period 24 hours. A more sophisticated
analysis or longer project could include collecting data from multiple buoys,
extended years, and/or additional information on storms
https://www.ncdc.noaa.gov/stormevents/.

Figure 5:  A smoothed periodogram of Wind Speed for buoy $\\#46025$ off the
coast of Santa Monica. The dashed lines are at 1/12 and 1/24, representing
12-hour and 24-hour periodicity.

###  3.5 Other Dynamic Datasets

There are myriad sources of dynamic data which are public and accessible. We
list a few additional resources here.

####  3.5.1 Baseball Data

Many students are interested in sports statistics. In particular,
statisticians and sabermatricians have worked over the years to compile
datasets and to think carefully about statistical methods applied to baseball
data. Today’s datasets can answer questions like: What is the probability of a
particular event happening given certain real time situations?

Albert (2010) provides curated seasonal batting data from 1871 to 2009 at
http://www.amstat.org/publications/jse/v18n3/mlb_batting.dat.txt. With the
passage of time, the examples of Derek Jeter and Alex Rodriguez as current
players and batting trends up to 2009 are less relevant to students interested
in sports, and what was once a current research question has become stagnant
and ossified. However, Albert (2010) provides sufficient detail about the
collection of files on the Lahman Baseball Database
(http://www.seanlahman.com/baseball-archive/statistics/ with data available as
a .zip file which can be read directly into R
http://seanlahman.com/files/database/lahman-csv_2015-01-24.zip) and tasks
required to filter and merge the master and batting files that the research
questions and teaching notes can be updated to reflect current players and
trends. Also see Jim Albert’s blog Exploring Baseball Data with R,
https://baseballwithr.wordpress.com/.

####  3.5.2 Cherry Blossom Ten Mile Run Data

Each year in April, the Cherry Blossom Ten Mile Run happens in Washington, DC
(http://www.cherryblossom.org/theraces/tenmile.php). The race results
(including finishing time, name, age, hometown, and pace) are available from
1999 to the present (although admittedly, the results after 2012 are more
difficult to scrape off the web due to a change in the format of how the
results are now posted). Kaplan and Nolan (2015) provide step-by-step
instructions (http://rdatasciencecases.org/) for scraping the data and
wrangling it into a format which can be used to investigate race results over
time, by age, or by gender.

####  3.5.3 Fatal Accidents and US Census Data

The National Highway Traffic Safety Administration collects data on all
fatalities suffered in motor vehicle traffic crashes. The data is posted
yearly and publicly as part of the Fatal Accident Reporting System,
http://www.nhtsa.gov/FARS. Each year is posted as a separate set of files,
available in different formats. For example, download the SAS data (readable
into R via the read.sas7bdat function in the `sas7bdat` package, the
sasxport.get function in the `Hmisc` package, or read.SAS in the `haven`
package) at ftp://ftp.nhtsa.dot.gov/fars/2014/SAS/.

Another dataset which is easily accessible is the census data
(http://factfinder2.census.gov/faces/nav/jsf/pages/index.xhtml). By looking at
the Decennial Census data and the Profile of General Population and Housing
Characteristics, students can study the ages and other demographic
characteristics of the general population as compared to those individuals
involved in fatal accidents333Original idea for this example from Laura
Kapitula, Grand Valley State University..

####  3.5.4 Climate Data

Witt (2013) describes using data in class to reveal important insights on
climate for statistics students. The data and analyses describe both decline
of Arctic sea ice and global temperature increase. The (static) data used in
the analyses is available through The Journal of Statistics Education,
http://www.amstat.org/publications/jse/v21n1/witt.pdf, but the authors also
provide links and information about the original data sources which are
dynamic.

Notably, realclimate.org provides a catalogue of many different types of
climate data and relevant source information
(http://www.realclimate.org/index.php/data-sources/). Witt (2013) also
reports:

> NASA’s Global Change Master Directory is available at
> http://gcmd.gsfc.nasa.gov/. NOAA’s National Climate Data Center maintains an
> extensive data directory available at http://www.ncdc.noaa.gov/oa/ncdc.html.
> Yet another good climate data source is the Data Guide maintained by the
> National Center for Atmospheric Research Climate at
> https://climatedataguide.ucar.edu/.

####  3.5.5 Iowa Liquor Sales Data

Iowa recently released an 800MB+ dataset containing all the weekly liquor
sales from January 1, 2014 to the present (https://data.iowa.gov/Economy/Iowa-
Liquor-Sales/m3tr-qhgy). The data seems to be updated monthly. Dan Nguyen has
provided some scripting code and SQL analysis for downloading and wrangling
the data (https://gist.github.com/dannguyen/18ed71d3451d147af414). Because the
data are available in a csv format, they can be read into R directly using
read.csv from the data.iowa.gov URL. The SQL code from the Nguyen’s github
site can be translated directly into R using the `dplyr` package. Note that
the rectangular structure and comma delimited format of the data make it
tempting to work with. The challenge for this particular dynamic dataset comes
with its size, which can make it unwieldy on many computers.

####  3.5.6 Medicare Inpatient Charges Data

Medicare inpatient charge data is available at https://www.cms.gov/Research-
Statistics-Data-and-Systems/Statistics-Trends-and-Reports/Medicare-Provider-
Charge-Data/Inpatient2013.html which links directly to a zipped csv file. An
analysis of the costs for Medicare Severity Diagnosis Related Group (MS-DRG)
(i.e., the procedure) can be done over multiple different covariates.

####  3.5.7 World Bank Data

The R package wbstats connects directly to the World Bank API, providing a
trove of global economic data. The package documentation and vignette provide
an easy start for students to download basic data of their choosing,
https://cran.r-project.org/web/packages/wbstats/vignettes/Using_the_wbstats_package.html.

####  3.5.8 Energy Consumption Data

The US Department of Energy collates energy usage at http://www.eia.gov/. For
example, at http://www.eia.gov/coal/data.cfm they provide complete data on
“coal production, consumption, exports, imports, stocks, mining, and prices."

####  3.5.9 Government Survey Data

Anthony Damico has compiled and documented a repository of data from many
dozens of different government surveys http://www.asdfree.com/, including
NHANES. Additionally, Damico espouses the merits of reproducible analysis done
in R with GitHub. As one example, the General Social Survey (GSS) contains
information on what Americans think of policies, issues, and priorities in the
US (http://gss.norc.org/).

###  3.6 Using the Examples

We try to communicate to the students that there is information in most data
sources. We want to be wary and attentive to issues of experimental design and
systematic biases. However, we do not want to leave our students feeling stuck
every time they encounter a dataset which has not been gathered from a large
randomized trial. Instead, we try to think of the pieces of the analysis that
can elicit information which is interesting or possibly hypothesis generating.

Students at every level are ready to examine complicated relationships using
tools which are accessible before taking multiple advanced statistics courses.
“Unfortunately, many instructors teach the sections on data analysis as
descriptive statistics, perhaps because this is what they experienced in their
first course. They emphasize the process of calculating numerical summaries
and making graphical displays, rather than using these as tools to explore
what the data are saying." (Notz, 2015)

To keep our students from getting stuck, “we must ask questions such as
whether the data allow generalization to a larger population, whether their
structure can be meaningfully described with the models we wish to fit, and
whether important subgroups or individuals were excluded from the data.
Exceptions, anomalies, outliers, and subgroups are best recognized and
understood in the context of the question being addressed." (De Veaux and
Velleman, 2015)

##  4  Conclusion

The few examples provided here give a glimpse into how to incorporate real and
dynamic data into introductory statistics and and courses beyond. The
Guidelines for Assessment and Instruction in Statistics Education (GAISE)
College Report implores instructors to make one small change in our courses to
keep up with changing technology and data (Carver _et al._ , 2016). In this
manuscript we have provided tools for helping educators seamlessly incorporate
dynamic data into a standard and modern introductory (or beyond) statistics
course. We also hope that statistics students will see the different ways of
accessing data and feel empowered to gather more information on their own. Our
courses must remain relevant to both the student experience and the research
community or we run the risk of become irrelevant. Small steps will get us
where we need to be.

##  5 Acknowledgements

This work was supported by Project MOSAIC (NSF grant 0920350, Phase II:
Building a Community around Modeling, Statistics, Computation, and Calculus).
Additionally, the data sources were found with help from many individuals. For
example, see Wes Stevenson’s great blog on importing data into R
(http://statistical-research.com/importing-data-into-r-from-different-
sources/). An important step forward for accessing data is Jenny Bryan’s R
package `googlesheets` (http://blog.revolutionanalytics.com/2015/09/using-the-
googlesheets-package-to-work-with-google-sheets.html). Many of the data sets
were inspired by work by or conversations with Nick Horton, Ben Baumer, Gabe
Chandler, Scott Grimshaw, Laura Kapitula, Danny Kaplan, Randy Pruim, Maddi
Cowen, Ciaran Evans, Samantha Morrison, and Janie Neal.

##  6 Appendix

The Appendix contains technical details of interest to students who want to
learn more about directly downloading data. R Studio allows for the teacher to
provide template R Markdown files which download data directly from the
Internet into the R Studio environment. Additionally, reproducible R Markdown
files teach good science and analysis. Horton _et al._ (2015) discuss the
merits of reproducible research as well as the ease of using R Markdown in the
classroom.

As part of the supplementary materials
(https://github.com/hardin47/DynamicData), we have provided R Markdown files
giving both R code and related pedagogical commentary associated with
different dynamic data examples for use in a statistics classroom. The data
have all been collected directly from outside sources - websites that are
updated periodically (and with information the students can access on their
own).

###  6.1 What is an API?

An Application Programming Interface (API) is a set of programming
instructions (written in code, giving the appropriate algorithm) for
downloading data from a website that contains - typically - vast amounts of
data. For example, Twitter has an API (https://dev.twitter.com/overview/api)
which tells programmers how to access tweets (and related information)
directly from the Twitter website. At any level, but particularly useful at an
introductory level, it is recommended to use an R (or other software) package
or function that allows R to speak to the API in order to download the data.
For example, the R package `wbstats` connects directly from R to the World
Bank API and `twitteR` provides an R interface to the Twitter API. Though the
examples provided do not generally rely on APIs or a related R interface, it
is good to be aware of APIs (and to tell your students!) in order to greatly
increase the sources of data available to you and your students.

###  6.2 R Markdown

The R Markdown files provided as supplementary materials make it easy for
instructors and students to perform reproducible analyses on data collection
through analysis and synthesis of the results. The students get the important
practice of tracking each and every step of their work. Reproducible analysis
has recently gotten a lot of press (Johnson, 2014; Stodden _et al._ , 2014),
and R Studio (RStudio Team, 2015) has produced a user friendly format for
combining R code with html (or LaTeX) word processing that has been used with
introductory statistics (Baumer _et al._ , 2014).

R Markdown has a short learning curve and will be straightforward for your
students to implement. One of the key considerations is that when running an R
Markdown file, the file does not pay attention to anything running locally. R
Markdown essentially restarts R, so the current state of your R session is
totally irrelevant. Recall that the purpose of R Markdown is to create
reproducible files, so the Markdown file should run on any computer anywhere
(regardless of the local environment).

Important note: Any package used in R Markdown needs to be installed prior to
compiling the .Rmd file. That is, for a line of code such as `library(dplyr)`,
before running the Markdown file, `install.packages("dplyr")` must be run from
within the console one time (ever).

###  6.3 Useful R functions

Some important R functions (in italics) and `packages` (in `typeface`) that
will help you and your students navigate importing and using data in R
include:

  * •

`dplyr` package for data wrangling in general; cheat sheet at
https://www.rstudio.com/wp-content/uploads/2015/02/data-wrangling-
cheatsheet.pdf

  * •

The glimpse function in `dplyr`

  * •

`tidyr` for converting between wide and long formats and for the very useful
extract_numeric() function (or readr::parse_numeric())

  * •

`ggplot2` for faceted graphing (also, `ggvis`) (Wickham and Sievert, 2016);
cheat sheet at https://www.rstudio.com/wp-
content/uploads/2015/03/ggplot2-cheatsheet.pdf

  * •

`openintro` (or packages that come with the textbooks you use) which are great
for pulling up any dataset from the text and building on it in class (Diez _et
al._ , 2012)

  * •

`mosaic` for consistent syntax and helpful functions used in introductory
statistics (Pruim _et al._ , 2014); cheat sheet at
https://cran.r-project.org/web/packages/mosaic/vignettes/MinimalR.pdf

  * •

`googlesheets` for loading data directly from Google spreadsheets

  * •

`lubridate` if you ever need to work with any date fields

  * •

`stringr` for text parsing and manipulation

  * •

`rvest` for scraping data off the web; `readxl` for reading excel data

  * •

`readr` (and `fread` with data.table) for loading large datasets with default
stringsAsFactors = FALSE

  * •

`tables` for nice looking summary tables.

## References

  * Albert (2010) Albert J (2010).  “Baseball Data at Season, Play-by-Play, and Pitch-by-Pitch Levels.”  _Journal of Statistics Education_ , 18(3).
  * American Statistical Association Undergraduate Guidelines Workgroup (2014) American Statistical Association Undergraduate Guidelines Workgroup (2014).  “2014 curriculum guidelines for undergraduate programs in statistical science.”  _Technical report_ , American Statistical Association, Alexandria, VA.  URL http://www.amstat.org/education/curriculumguidelines.cfm.
  * Baumer _et al._ (2014) Baumer B, Çetinkaya-Rundel M, Bray A, Loi L, Horton N (2014).  “R Markdown: Integrating A Reproducible Analysis Tool into Introductory Statistics.”  _Technology Innovations in Statistics Education_ , 8(1).  URL http://escholarship.org/uc/item/90b2f5xh.
  * Brown and Kass (2009) Brown EN, Kass RE (2009).  “What is Statistics?”  _The American Statistician_ , 63(2), 105–110.
  * Carver _et al._ (2016) Carver R, Everson M, Gabrosek J, Rowell GH, Horton N, Lock R, Mocko M, Rossman A, Velleman P, Witmer J, Wood B (2016).  “Guidelines for Assessment and Instruction in Statistics Education (GAISE) College Report.”  _Technical report_ , American Statistical Association.  URL http://www.amstat.org/education/gaise/collegeupdate/GAISE2016_DRAFT.pdf.
  * Cobb (1991) Cobb GW (1991).  “Teaching Statistics: More Data, Less Lecturing.”  _UME Trends_ , pp. 3–7.
  * Cobb (1992) Cobb GW (1992).  “Teaching Statistics.”  _In Lynn A. Steen (ed), Heeding the call for change: suggestions for curricular action (MAA Notes No. 22)_ , pp. 3–43.
  * Cobb (2007) Cobb GW (2007).  “The Introductory Statistics Course: A Ptolemaic Curriculum?”  _Technology Innovations in Statistics Education_ , 1(1).  URL https://escholarship.org/uc/item/6hb3k0nz.
  * Cobb (2011) Cobb GW (2011).  “Teaching statistics: some important tensions.”  _Chilean Journal of Statistics_ , 2(1), 31–62.
  * De Veaux and Velleman (2015) De Veaux R, Velleman P (2015).  “Teaching Statistics Algorithmically or Stochastically Misses the Point: Why not Teach Holistically? (Online discussion of “Mere Renovation is Too Little Too Late: We Need to Rethink Our Undergraduate Curriculum From the Ground Up," by George W Cobb, The American Statistician).”  _The American Statistician_ , 69(4).
  * Diez _et al._ (2012) Diez DM, Barr CD, Çetinkaya-Rundel M (2012).  _openintro: OpenIntro data sets and supplemental functions_.  R package version 1.4, URL http://CRAN.R-project.org/package=openintro.
  * GAISE College Group (2005) GAISE College Group (2005).  “Guidelines for Assessment and Instruction in Statistics Education.”  _Technical report_ , American Statistical Association.  URL http://www.amstat.org/education/gaise.
  * Gould (2010) Gould R (2010).  “Statistics and the Modern Student.”  _International Statistical Review_ , 78(2), 297–315.
  * Gould and Çetinkaya-Rundel (2013) Gould R, Çetinkaya-Rundel M (2013).  “Teaching Statistical Thinking in the Data Deluge.”  In T Wassong, D Frischemeier, PR Fischer, R Hochmuth, P Bender (eds.), _Using Tools for Learning Statistics and Mathematics_ , pp. 377–391. Springer.
  * Grimshaw (2015) Grimshaw S (2015).  “A Framework for Infusing Authentic Data Experiences Within Statistics Courses.”  _The American Statistician_ , 69(4), 307–314.  URL http://arxiv.org/abs/1507.08934.
  * Horton _et al._ (2015) Horton NJ, Baumer B, Wickham H (2015).  “Setting the stage for data science: integration of data management skills in introductory and second courses in statistics.”  _CHANCE_ , 28(2), 40–50.  URL http://arxiv.org/abs/1502.00318.
  * Johnson (2014) Johnson G (2014).  “New Truths That Only One Can See.”  URL http://www.nytimes.com/2014/01/21/science/new-truths-that-only-one-can-see.html.
  * Kaplan and Nolan (2015) Kaplan D, Nolan D (2015).  “Modeling Runners’ Times in the Cherry Blossom Race.”  In _Data Science in R: A Case Studies Approach to Computational Reasoning and Problem Solving_ , pp. 45–104. CRC Press.
  * Kuiper and Sturdivant (2015) Kuiper S, Sturdivant RX (2015).  “Using Online Game-Based Simulations to Strengthen Students’ Understanding of Practical Statistical Issues in Real-World Data Analysis.”  _The American Statistician_ , 69(4), 354–361.
  * Notz (2015) Notz W (2015).  “Vision or Bad Dream? (Online discussion of “Mere Renovation is Too Little Too Late: We Need to Rethink Our Undergraduate Curriculum From the Ground Up,” by George W Cobb, The American Statistician).”  _The American Statistician_ , 69(4).
  * Pruim _et al._ (2014) Pruim R, Kaplan D, Horton N (2014).  _mosaic: Project MOSAIC (mosaic-web.org) Statistics and Mathematics Teaching Utilities_.  R package version 0.9-1-3, URL http://CRAN.R-project.org/package=mosaic.
  * RStudio Team (2015) RStudio Team (2015).  _RStudio: Integrated Development Environment for R_.  RStudio, Inc., Boston, MA.  URL http://www.rstudio.com/.
  * Stodden _et al._ (2014) Stodden V, Leisch F, Peng RD (eds.) (2014).  _Implementing Reproducible Research_.  Chapman and Hall / CRC Press.
  * Wickham (2014) Wickham H (2014).  “Tidy Data.”  _Journal of Statistical Software_ , 59(10).  URL http://www.jstatsoft.org/v59/i10/.
  * Wickham and Sievert (2016) Wickham H, Sievert C (2016).  _ggplot2: Elegant Graphics for Data Analysis_.  Springer New York.  URL http://had.co.nz/ggplot2/book.
  * Witt (2013) Witt G (2013).  “Using Data from Climate Science to Teach Introductory Statistics.”  _Journal of Statistics Education_ , 21(1).
  * Workgroup on Undergraduate Statistics (2000) Workgroup on Undergraduate Statistics (2000).  “Guidelines for Undergraduate Statistics Programs, http://www.amstat.org/education/curriculumguidelines.cfm, accessed August 18, 2013.”  _Technical report_ , American Statistical Association.
  * Zhu _et al._ (2013) Zhu Y, Hernandez LM, Mueller P, Dong Y, Forman MR (2013).  “Data Acquisition and Preprocessing in Studies on Humans: What is Not Taught in Statistics Classes?”  _The American Statistician_ , 67, 235–241.


# Dynamical correlation energy of metals in large basis sets from downfolding
and composite approaches

James M. Callahan  Department of Chemistry, Columbia University, New York, New
York 10027 USA  Malte F. Lange  Department of Chemistry, Columbia University,
New York, New York 10027 USA  Timothy C. Berkelbach  Department of Chemistry,
Columbia University, New York, New York 10027 USA  Center for Computational
Quantum Physics, Flatiron Institute, New York, New York 10010 USA
<EMAIL_ADDRESS>
###### Abstract

Coupled-cluster theory with single and double excitations (CCSD) is a
promising ab initio method for the electronic structure of three-dimensional
metals, for which second-order perturbation theory (MP2) diverges in the
thermodynamic limit. However, due to the high cost and poor convergence of
CCSD with respect to basis size, applying CCSD to periodic systems often leads
to large basis set errors. In a common “composite” method, MP2 is used to
recover the missing dynamical correlation energy through a focal-point
correction, but the inadequacy of MP2 for metals raises questions about this
approach. Here we describe how high-energy excitations treated by MP2 can be
“downfolded” into a low-energy active space to be treated by CCSD. Comparing
how the composite and downfolding approaches perform for the uniform electron
gas, we find that the latter converges more quickly with respect to the basis
set size. Nonetheless, the composite approach is surprisingly accurate because
it removes the problematic MP2 treatment of double excitations near the Fermi
surface. Using the method to estimate the CCSD correlation energy in the
combined complete basis set and thermodynamic limits, we find CCSD recovers
over 90% of the exact correlation energy at $r_{s}=4$. We also test the
composite and downfolding approaches with the random-phase approximation used
in place of MP2, yielding a method that is more effective but more expensive.

## Introduction

Ground-state electronic properties of metallic solids have traditionally been
computed using density functional theory (DFT) [1, 2, 3, 4], which is
partially justified by the fact that many popular functionals are
parameterized by numerically exact results on the uniform electron gas (UEG)
[5, 6, 7, 8]. In recent years, interest has grown around the application of ab
initio, wavefunction-based electronic structure techniques for condensed-phase
systems [9, 10, 11, 12, 13, 14, 15, 16], since they do not suffer from
uncontrolled errors inherent to the DFT exchange-correlation functional [6,
17, 18]. Promising methods in this direction include the random-phase
approximation (RPA) [19, 20, 16] and coupled-cluster theory [21, 22, 23, 24,
25, 26]. Importantly, both of these methods preclude the well-known
divergences of finite-order perturbation theories, such as second-order
Møllet-Plesset perturbation theory (MP2), via an infinite-order resummation
[27, 28, 29, 30, 31].

Although coupled-cluster theory has been successfully applied to an increasing
number of atomistic semiconductors and insulators [32, 14, 15, 33, 34, 35, 36,
37], its applicability for metals has been primarily focused around the UEG,
also known as jellium [29, 38, 39, 40, 41, 42, 43, 44, 45]. Despite their
reasonable accuracy, these calculations have demonstrated the typical slow
convergence of the correlation energy as a function of the number of virtual
(unoccupied) orbitals included [11, 46, 47, 40]. This slow convergence is
especially problematic because of the high cost of coupled-cluster
calculations with large basis sets. For example, coupled-cluster theory with
single and double excitations (CCSD) has a computational cost that scales as
$O(N^{2}M^{4})$, where $N$ and $M$ are the number of electrons and basis
functions respectively. To date, results near the complete basis set (CBS)
limit have been primarily computed via the extrapolation of results obtained
with a finite, increasing number of basis functions [40, 45, 24], although
explicitly correlated [48] and transcorrelated [49, 50] methods provide
promising alternative approaches.

Composite methods (sometimes called focal point methods) are a simple,
alternative class of approaches for recovering dynamical correlation within
large basis sets [51, 52, 53]. A common composite scheme combines the results
of high-level and low-level theories using three calculations. For example,
using CCSD as the high-level theory and MP2 as a low-level theory, the CCSD
correlation energy in a large basis is approximated as

$E_{\mathrm{CCSD}}(M)\approx
E_{\mathrm{CCSD}}(M_{\mathrm{act}})+E_{\mathrm{MP2}}(M)-E_{\mathrm{MP2}}(M_{\mathrm{act}})$
(1)

where $M_{\mathrm{act}}<M$ is the number of “active” basis functions. Refs.
54, 55 provide a similar but more sophisticated CCSD/MP2 composite method,
based on an analysis of the basis set convergence of various diagrammatic
contributions to the correlation energy. While such CCSD/MP2 composite
approaches have been applied successfully to a number of semiconductors and
insulators [54, 56, 36, 57], their applicability to metals is questionable
because of the failures of MP2 theory. One goal of this work is to test the
composite CCSD/MP2 approach for metals.

A more theoretically satisfying approach would be to perform a single
calculation where low-energy excitations near the Fermi surface are treated
with CCSD and are coupled to high-energy excitations treated with MP2. This
particular approach, which is similar to tailored CC [58] and the broader
class of active-space CC methods [59, 60, 61, 62, 63], has variously been
called CC/PT [64], CCSD-MP2 [65], and multilevel CC [66, 67]. Two of us
(M.F.L. and T.C.B.) recently tested this method for a few simple atomistic
semiconductors and insulators [56], and here we aim to assess its performance
for metals, where the differences between CCSD and MP2 are more striking.
Since the effects of the frozen high-energy MP2 amplitudes are folded down
onto the low-energy CCSD amplitudes (see below), we refer to this method as a
“downfolding” approach. In principle, this downfolding CCSD/MP2 method should
provide a distinct advantage over the conceptually simpler composite approach,
as downfolding does not include the MP2 treatment of low-energy excitations
that are responsible for divergence in the thermodynamic limit (TDL). After
providing theoretical details of these two methods, we compare their
performance for the UEG at a fixed number of electrons and in the TDL. Before
concluding, we also examine the straightforward use of the RPA in place of
MP2.

## Theory

Here we briefly review the theory underlying the downfolding and composite
approaches. The $N$ occupied spin-orbitals are indexed by $i,j,k,l$; the
$\left(M-N\right)$ virtual orbitals by $a,b,c,d$; and the $M$ general orbitals
by $p,q,r,s$. The MP2 and coupled-cluster with double excitations (CCD)
correlation energies are given by

$E_{\mathrm{c}}=\frac{1}{4}\sum_{ijab}t_{ij}^{ab}\langle ij||ab\rangle$ (2)

where $t_{ij}^{ab}$ are amplitudes of the double excitation operator
$T_{2}=\frac{1}{4}\sum_{ijab}t_{ij}^{ab}a_{a}^{\dagger}a_{b}^{\dagger}a_{j}a_{i}$
and $\langle pq||rs\rangle$ are antisymmetrized two-electron repulsion
integrals; contributions from single excitations vanish because the UEG has no
capacity for orbital relaxation by symmetry. At lowest order in perturbation
theory,

$t_{ij}^{ab}=\frac{\langle
ab||ij\rangle}{\varepsilon_{i}+\varepsilon_{j}-\varepsilon_{a}-\varepsilon_{b}},$
(3)

and Eq. (2) gives the MP2 correlation energy. The high density of states at
the Fermi surface and the long-ranged nature of the Coulomb potential are
together responsible for the divergence of MP2 in the TDL. By contrast, the
CCD amplitudes solve a system of nonlinear equations

$0=\langle\Phi_{ij}^{ab}|e^{-T_{2}}He^{T_{2}}|\Phi\rangle,$ (4)

where $H$ is the electronic Hamiltonian. A standard approach for reaching the
CBS limit is to perform a series of calculations with increasing $M$ and use a
$M^{-1}$ extrapolation.

In both the composite and downfolding approaches, we partition the orbitals
into a set of $M_{\mathrm{act}}$ active orbitals, composed of all occupied
orbitals and the low-energy virtual orbitals, and a set of
$\left(M-M_{\mathrm{act}}\right)$ frozen (inactive) orbitals, composed of the
high-energy virtual orbitals. In principle, occupied orbitals can also be
partitioned, but typically they do not significantly contribute to the
computational cost. In the composite CCD/MP2 approach, the correlation energy
is calculated according to Eq. (1). Importantly for metals, the low-energy
active space double excitations are treated by CCD and not by MP2, so we
expect the method to be well-behaved in the thermodynamic limit.

In the downfolding CCD/MP2 approach, the double excitation operator $T_{2}$ is
partitioned into internal excitations fully contained within the active space
and external excitations that involve at least one frozen orbital,
$T_{2}=T_{2}^{\mathrm{(int)}}+T_{2}^{\mathrm{(ext)}}$. Fixing the
$T_{2}^{\mathrm{(ext)}}$ amplitudes to their MP2 values via Eq. (3), the
downfolding method involves first solving Eqs. (4) for only the internal
amplitudes and then evaluating the correlation energy expression Eq. (2) using
both the internal and external amplitudes. Compared to the $O(N^{2}M^{4})$
cost of full CCD, the composite approach has
$O(N^{2}M_{\mathrm{act}}^{4})+O(N^{2}M^{2})$ cost and the downfolding approach
has $O(N^{2}M_{\mathrm{act}}^{2}M^{2})$ cost, which can provide significant
savings, depending on the practical value of the ratio $M_{\mathrm{act}}/M$.

Let us now provide more insight into the “downfolding” perspective. Note that,
because the internal and external excitation operators commute, the defining
energy and amplitude equations of the downfolding approach can also be written

$\displaystyle\begin{split}E_{\mathrm{c}}&=\langle\Phi|e^{-T_{2}^{\mathrm{(int)}}}(\bar{H}-E_{\mathrm{HF}})e^{T_{2}^{\mathrm{(int)}}}|\Phi\rangle\\\
&=E_{\mathrm{MP2}}^{\mathrm{(ext)}}+\frac{1}{4}\sum_{ijab}^{\mathrm{active}}t_{ij}^{ab}\langle
ij||ab\rangle\end{split}$ (5a) $\displaystyle 0$
$\displaystyle=\langle\Phi_{ij}^{ab}|e^{-T_{2}^{\mathrm{(int)}}}\bar{H}e^{T_{2}^{\mathrm{(int)}}}|\Phi\rangle\quad(i,j,a,b)\
\mathrm{active}$ (5b)

where $\bar{H}=e^{-T_{2}^{\mathrm{(ext)}}}He^{T_{2}^{\mathrm{(ext)}}}$ (with
fixed $T_{2}^{\mathrm{(ext)}}$ as detailed above) and
$E_{\mathrm{MP2}}^{\mathrm{(ext)}}=\frac{1}{4}\sum_{ijab}^{\mathrm{(ext)}}t_{ij}^{ab}\langle
ij||ab\rangle$ is the MP2 correlation energy due to external excitations.
These resemble ordinary CCD energy and amplitude equations within the active
space only, except that the bare Hamiltonian $H$ is replaced by an effective
Hamiltonian $\bar{H}$ that is similarity-transformed by the external
excitation amplitudes. The effective Hamiltonian within the active space can
be expressed as

$\begin{split}\bar{H}-E_{\mathrm{HF}}&=E_{\mathrm{MP2}}^{\mathrm{(ext)}}+\sum_{pq}^{\textrm{active}}F_{pq}\\{a_{p}^{\dagger}a_{q}\\}\\\
&\hskip
10.00002pt+\frac{1}{4}\sum_{pqrs}^{\textrm{active}}W_{pqrs}\\{a_{p}^{\dagger}a_{q}^{\dagger}a_{s}a_{r}\\}+\dots,\end{split}$
(6)

where $\\{\dots\\}$ indicates normal ordering of the operators. This effective
Hamiltonian can be seen to contain effective one- and two-body interactions
that are frequency independent [21], in contrast to other downfolding
approaches like the constrained random-phase approximation [68, 69]. For
example, the all-occupied two-body interaction becomes

$W_{ijkl}=\langle
ij||kl\rangle+\frac{1}{2}\sum_{ab}{\vphantom{\sum}}^{\prime}\frac{\langle
ij||ab\rangle\langle
ab||kl\rangle}{\varepsilon_{k}+\varepsilon_{l}-\varepsilon_{a}-\varepsilon_{b}},$
(7)

where the primed summation indicates that one or both of $a,b$ are inactive
virtual orbitals. The frequency independence can be understood because our
observable is the total energy rather than a spectral function. To summarize,
an approach that solves the internal CCD amplitude equations in the presence
of frozen external amplitudes is equivalent to a CCD calculation in an active
space of orbitals using an effective (downfolded) Hamiltonian that is
similarity-transformed by the external excitation operator.

It is straightforward to show that the composite approach, normally understood
as a three-step procedure as shown in Eq. (1), is equivalent to Eqs. (5) but
where the effective Hamiltonian $\bar{H}$ is replaced by the bare Hamiltonian
$H$ in the amplitude equations (5b). From this perspective, the performance
differences between the downfolding and composite approaches are attributable
to the screening of the integrals in the effective Hamiltonian when
determining the internal amplitudes. Nevertheless, we reiterate that the
composite CCD/MP2 approach is expected to perform well because it replaces the
problematic MP2 treatment of low-energy, internal double excitations with a
well-behaved CCD treatment.

Figure 1: Basis set convergence of the correlation energy of the $r_{s}=4$ UEG
with $N=332$ electrons for MP2 (squares), CCD (circles), composite CCD/MP2
(crosses), and downfolding CCD/MP2 (diamonds). The top axis shows the
percentage of virtual orbitals that are active for the composite and
downfolding methods, compared to the “target” calculation with $M=1502$. Both
the composite and downfolding methods interpolate between the “target” MP2
calculation (leftmost square) and the “target” CCD calculation (leftmost
circle).

## Results and discussion

We study the UEG as the simplest model of metals. A brief review of the UEG
model in finite cells with finite plane-wave basis sets is given in the
Appendix and we refer the reader to the literature for more details [31, 70,
40]. To illustrate the performance of the composite and downfolding methods,
we focus on the Wigner-Seitz radius $r_{s}=4$ (corresponding to the
approximate valence electron density of metallic sodium), where CCD has been
found to recover about 85% of the correlation energy [40, 43]. We use a
twisted boundary condition by performing calculations at the Baldereschi point
[71], which has been shown to provide smoother convergence to the TDL [72,
73].

In Fig. 1, we show basis set convergence of the correlation energy for a
finite UEG with $N=332$ electrons. The uncorrected MP2 and CCD correlation
energies exhibit their typical slow convergence with increasing basis set size
and show asymptotic behavior where the basis set error decays as $M^{-1}$.
Extrapolation to the CBS limit yields
$E_{\mathrm{c}}/N=-0.0401~{}E_{\mathrm{h}}$ for MP2 and
$E_{\mathrm{c}}/N=-0.0262~{}E_{\mathrm{h}}$ for CCD. At the largest finite
basis shown, $M=1502$, the results exhibit a significant basis set error of
about $0.01~{}E_{\mathrm{h}}$ for both methods, highlighting the challenge of
recovering dynamical correlation in metals with large basis sets. Importantly,
we emphasize that the MP2 correlation energy does not diverge for any finite
system but only upon extrapolation to the TDL (see below).

Recall that the CCD/MP2 composite and downfolding approaches involve both a
“target” number of orbitals $M$ and an active number of orbitals
$M_{\mathrm{act}}$. In Fig. 1, we show results obtained for $M=1502$ as
$M_{\mathrm{act}}$ is varied. By construction, both methods yield the target
MP2 correlation energy when there are no active virtual orbitals
($M_{\mathrm{act}}=N$) and the target CCD correlation energy when all orbitals
are active ($M_{\mathrm{act}}=M$). We observe that both methods converge
smoothly to the target CCD result and that the downfolding approach exhibits a
faster convergence, due to its coupling between the internal and external
excitation spaces. We also see similar behavior for other numbers of electrons
and densities (not shown), indicating that neither finite-size effects nor the
specific metallic density changes the overall picture.

Figure 2: Absolute error in the correlation energy for the data in Figure 1,
shown on a logarithmic scale, relative to the “target” CCD result with
$M=1502$. The top axis and symbols have the same meaning as Figure 1. Dotted
black lines are shown as a guide for various power law exponents $\alpha$ as
discussed in the text.

To better quantify the rate of convergence, in Fig. 2 we plot the absolute
deviation of the correlation energy from the “target” CCD result obtained with
$M=1502$. The error is plotted as a function of the difference between the
inverse number of active orbitals and the inverse number of total orbitals and
analyzed in terms of the power law $|\Delta
E_{\mathrm{c}}|\propto\left[M_{\mathrm{act}}^{-1}-M^{-1}\right]^{\alpha}$. We
compare the convergence of traditional CCD, the composite approach, and the
downfolding approach. For plain CCD, we see linear convergence of the
correlation energy, with $\alpha\approx 1$, over a large range of basis set
sizes. The composite method exhibits an early, rapid convergence reaching a
maximum scaling of around $\alpha\approx 2$ before slowing to the same
$\alpha\approx 1$ convergence as $M_{\mathrm{act}}$ approaches $M$. The rapid
convergence for small $M_{\mathrm{act}}$ is responsible for absolute errors
that are about one order of magnitude better than those obtained by simple
truncation. In fact, the plain CCD result does not obtain m$E_{\mathrm{h}}$
accuracy until essentially all orbitals are correlated, whereas the composite
result achieves this accuracy when only 50% of the virtual orbitals are
correlated in the expensive CCD calculation; this results in a speedup of a
factor of 16 compared to the full CCD calculation. Finally, the downfolding
result exhibits rapid but non-monotonic convergence, making it difficult to
extract a power law. Before slightly overshooting the “target”, the power law
exponent reaches $\alpha\approx 3$ or better, a significant improvement over
the composite CCD/MP2 and standard CCD approaches. This fast rate of
convergence provides m$E_{\mathrm{h}}$ accuracy when about one third of the
virtual orbitals are correlated, giving a speedup of a factor of 9.

Figure 3: Thermodynamic limit convergence of the correlation energy of the
$r_{s}=4$ UEG for systems with $N=90$–2392 electrons. (a) MP2 results for
basis set sizes indicated (squares) and in the extrapolated CBS limit (stars).
Thin diamonds show the CBS MP2 correction, $\delta^{(2)}\left(M/N\approx
4.0\right)$, indicated by the double-headed arrow at $N=210$. (b) CCD
(circles) and composite CCD/MP2 (crosses) results at the same (active) basis
set sizes. For composite CCD/MP2, we applied the CBS MP2 correction
$\delta^{(2)}\left(M\right)$. The grey line shows our TDL extrapolation using
the largest five systems for $M/N\approx 4.0$.

The good performance of the composite approach indicates that MP2, while an
inapplicable theory for three-dimensional metals, is safe to use for basis set
corrections. As discussed above, the reason for this applicability can be
understood by considering the MP2 correction that is applied in Eq. (1). This
correction is a difference between two MP2 correlation energies, both of which
correlate orbitals near the Fermi surface. These two MP2 energies are
separately divergent in the TDL, but their difference is not; moreover, this
difference is precisely $E_{\mathrm{MP2}}^{(\mathrm{ext})}$ defined
previously.

The reliable scaling of the MP2 basis set error at large $M$ suggests that the
MP2 CBS limit can be obtained by extrapolation for any given number of
electrons. In contrast, the asymptotic scaling regime for the CCD correlation
energy cannot always be reached. Thus, for any calculation performed with a
given $M$, we propose to add the MP2 correlation energy difference
$\delta^{(2)}\left(M\right)=E_{\mathrm{MP2}}\left(\infty\right)-E_{\mathrm{MP2}}\left(M\right)$,
where $E_{\mathrm{MP2}}\left(\infty\right)$ is obtained by $M^{-1}$
extrapolation. Having obtained an estimate of the CBS limit for a given number
of electrons, the finite-size extrapolation to the TDL can be done separately.
We expect this scheme to be not only more reliable than extrapolating CCD on
its own but also less costly since it involves more MP2 and fewer CCD
calculations.

In Fig. 3(a), we plot the MP2 correlation energy as a function of the inverse
number of electrons for finite UEG systems containing $N=90$–2392 electrons.
For each system size, we performed MP2 calculations at different basis set
sizes (grey squares), and the top four grey curves connect systems with a
similar ratio of $M/N$. We then performed $M^{-1}$ extrapolations at each
particle number using data from $M/N\approx 3.0,3.5,4.0,4.5$ to obtain the CBS
limit at each system size (black stars). On approach to the TDL, the MP2
correlation energy diverges, as seen most easily in our largest calculations
for the smaller values of $M/N$. Despite the divergence of its components, the
MP2 CBS correction $\delta^{(2)}$, plotted explicitly in Fig. 3(a) for
$M/N\approx 4$ (pink thin diamonds), does not diverge and is thus safe for use
in metallic systems.

To confirm that the success of the composite CCD/MP2 method holds on approach
to the TDL, in Fig. 3(b) we plot results for CCD (green circles) and composite
CCD/MP2 (orange crosses). For the composite CCD/MP2 result, we applied the MP2
CBS correction $\delta^{(2)}\left(M\right)$, although we note that a similar
CBS correction could also be combined with the downfolding CCD/MP2 approach.
In contrast to MP2, CCD is well-defined in the TDL; however, at each system
size, convergence to the CBS limit is slow, as shown by the green data sets.
In contrast, we observe that composite CCD/MP2 has much faster convergence to
the TDL, as shown by the orange data sets. Performing a $N^{-1}$ extrapolation
of the CBS-corrected $M/N\approx 4.0$ results gives us a CBS and TDL
extrapolated correlation energy of $E_{\mathrm{c}}=-0.0293~{}E_{\mathrm{h}}$,
in general agreement with past CCD results[40]. For comparison, the exact
value [74] is $E_{\mathrm{c}}=-0.0318~{}E_{\mathrm{h}}$, indicating that CCD
recovers over 90% of the correlation energy. The same finite-size
extrapolation of the CBS-corrected $M/N\approx 3.0$ data gives a correlation
energy that differs by only $0.7~{}\mathrm{m}E_{\mathrm{h}}$, indicating the
excellent convergence of the composite method.

Figure 4: The same as in Fig. 1, except for $r_{s}=1$. In addition to the
methods shown in Fig. 1, we also include RPA (pentagons), composite CCD/RPA
(thin crosses), and downfolding CCD/RPA (pluses).

Before concluding, we recognize that a variety of other low-level theories can
be combined with CCD, in both the composite and downfolding manner. Of
particular interest is the RPA, which is more appropriate for metals than MP2
but also more expensive. The so-called direct RPA is most promising
computationally due to its low cost when two-electron integrals are handled by
density fitting. Here we instead test the full particle-hole RPA, where the
$T_{2}$ amplitudes maintain proper anti-symmetry so that the theory is free of
self-interaction error. The RPA amplitudes are the solution of the CCD
equations where only selected terms are retained [75, 70]. In Fig. 4, we show
the same results as in Fig. 1, except for the electron density corresponding
to $r_{s}=1$ (due to RPA convergence problems at larger values of $r_{s}$
[70]). We observe the same overall trends as we did at $r_{s}=4$ for the MP2
and CCD calculations with finite basis sets. However, the RPA results follow
the CCD results much more closely, which also makes composite and downfolding
CCD/RPA significantly outperform the corresponding MP2 methods.
Unsurprisingly, the improved performance of the downfolding approach compared
to the composite approach is even more marginal than for CCD/MP2.

## Conclusions and outlook

We have described and analyzed two approaches for eliminating basis set error
in the CCD correlation energy of metals, using the simple UEG model. Our
results indicate that these methods allow for aggressive freezing of virtual
orbitals or approximation of external amplitudes, leading to significant
reductions in computational cost. Although the downfolding CCD/MP2 approach is
slightly more accurate, we find that the simpler composite CCD/MP2 approach is
surprisingly effective because divergent contributions near the Fermi surface
do not contribute to the basis set correction.

In this work, we have addressed post-Hartree-Fock basis set errors in the
canonical orbital basis, but the methods we presented could also be
straightforwardly applied in a basis of localized orbitals. For example, Refs.
76, 57 used CCSD/MP2 and CCSD/RPA composite approaches to mitigate basis set
errors in quantum embedding calculations. Localized orbital basis sets mix
orbitals near and far from the Fermi surface, so it will be interesting to
test how the composite and downfolding approaches perform when using these
localized orbitals for metals.

Future work will focus on applying these techniques to atomistic metals using
natural orbitals [77], to the excited-state properties of metals [41, 44], and
to higher-level theories of correlation [43]. For example, we imagine that a
composite CCSDT/CCSD or CCSDT/CCSD(T) approach would provide quantitative
accuracy for metals while precluding the failure [31] of the otherwise
successful treatment of perturbative triple excitations.

###### Acknowledgements.

We thank Verena Neufeld and Sandeep Sharma for comments on the manuscript.
This work was supported in part by the National Science Foundation Graduate
Research Fellowship under Grant No. DGE-1644869 (M.F.L.), by the Department of
Defense through the National Defense Science & Engineering Graduate (NDESG)
Fellowship Program (J.M.C.), and by the National Science Foundation under
Grant No. CHE-1848369 (T.C.B.). We acknowledge computing resources from
Columbia University’s Shared Research Computing Facility project, which is
supported by NIH Research Facility Improvement Grant 1G20RR030893-01, and
associated funds from the New York State Empire State Development, Division of
Science Technology and Innovation (NYSTAR) Contract C090171, both awarded
April 15, 2010. The Flatiron Institute is a division of the Simons Foundation.

## Appendix: Uniform electron gas

Working in a single-particle basis of plane waves with momenta
${\bm{k}}=(2\pi/L)(n_{x},n_{y},n_{z})$ and cell of volume $L^{3}$ with
periodic boundary conditions, the UEG Hamiltonian is given by

$\begin{split}H&=\sum_{{\bm{k}}\sigma}\frac{k^{2}}{2}a_{{\bm{k}}\sigma}^{\dagger}a_{{\bm{k}}\sigma}\\\
&+\frac{1}{2}\sum_{{\bm{k}}_{1}{\bm{k}}_{2}{\bm{k}}_{3}{\bm{k}}_{4}}^{\prime}\sum_{\sigma\sigma^{\prime}}\langle{\bm{k}}_{1}\sigma,{\bm{k}}_{2}\sigma^{\prime}|{\bm{k}}_{3}\sigma,{\bm{k}}_{4}\sigma^{\prime}\rangle
a^{\dagger}_{{\bm{k}}_{1}\sigma}a^{\dagger}_{{\bm{k}}_{2}\sigma^{\prime}}a_{{\bm{k}}_{4}\sigma^{\prime}}a_{{\bm{k}}_{3}\sigma}\end{split}$
(8)

where the primed summation requires
${\bm{k}}_{1}+{\bm{k}}_{2}={\bm{k}}_{3}+{\bm{k}}_{4}$. The two-electron
repulsion integrals are given by

$\langle{\bm{k}}_{1}\sigma,{\bm{k}}_{2}\sigma^{\prime}|{\bm{k}}_{3}\sigma,{\bm{k}}_{4}\sigma^{\prime}\rangle=v({\bm{k}}_{1}-{\bm{k}}_{3})\delta_{{\bm{k}}_{1}+{\bm{k}}_{2},{\bm{k}}_{3}+{\bm{k}}_{4}}$
(9)

where the Ewald potential is

$\displaystyle v({\bm{k}})$
$\displaystyle=\begin{cases}\dfrac{4\pi}{L^{3}k^{2}}&k\neq 0\\\ \\\
v_{\mathrm{M}}&k=0,\end{cases}$ (10)

and $v_{\mathrm{M}}=2.837297479/L$ is the Madelung constant of the cell [78].
The $N$-electron reference determinant has the lowest-energy $N/2$ plane wave
orbitals doubly occupied and the HF orbital energies are given by
$\varepsilon({\bm{k}})=k^{2}/2-v_{\mathrm{M}}\theta(k_{\mathrm{F}}-k)$, where
$k_{\mathrm{F}}$ is the Fermi momentum. We restrict our calculations to
closed-shell configurations, which allows only certain “magic numbers” of
electrons and orbitals. At the Baldereschi point, the first few magic numbers
are $2,8,14,22,34,40,52$.

## Data availability statement

The data that support the findings of this study are available from the
corresponding author upon reasonable request.

## References

  * Kohn and Sham [1965] W. Kohn and L. J. Sham, “Self-consistent equations including exchange and correlation effects,” Phys. Rev. 140, A1133–A1138 (1965).
  * Kresse and Furthmüller [1996] G. Kresse and J. Furthmüller, “Efficiency of ab-initio total energy calculations for metals and semiconductors using a plane-wave basis set,” Comput. Mater. Sci. 6, 15–50 (1996).
  * Burke [2012] K. Burke, “Perspective on density functional theory,” J. Chem. Phys. 136, 150901 (2012).
  * Jones [2015] R. O. Jones, “Density functional theory: Its origins, rise to prominence, and future,” Rev. Mod. Phys. 87, 897–923 (2015).
  * Ceperley and Alder [1980] D. M. Ceperley and B. J. Alder, “Ground state of the electron gas by a stochastic method,” Phys. Rev. Lett. 45, 566–569 (1980).
  * Huang and Carter [2008] P. Huang and E. A. Carter, “Advances in Correlated Electronic Structure Methods for Solids, Surfaces, and Nanostructures,” Annu. Rev. Phys. Chem. 59, 261–290 (2008).
  * Giuliani and Vignale [2005] G. Giuliani and G. Vignale,  _Quantum theory of the electron liquid_, 1st ed. (Cambridge University Press, 2005).
  * Loos and Gill [2016] P. F. Loos and P. M. Gill, “The uniform electron gas,” WIRES Comput. Mol. Sci. 6, 410–429 (2016).
  * Pisani  _et al._ [2008] C. Pisani, L. Maschio, S. Casassa, M. Halo, M. Schütz, and D. Usvyat, “Periodic local MP2 method for the study of electronic correlation in crystals: Theory and preliminary applications,” J. Comput. Chem. 29, 2113–2124 (2008).
  * Grüneis, Marsman, and Kresse [2010] A. Grüneis, M. Marsman, and G. Kresse, “Second-order Møller-Plesset perturbation theory applied to extended systems. II. Structural and energetic properties,” J. Chem. Phys. 133, 074107 (2010).
  * Shepherd  _et al._ [2012a] J. J. Shepherd, A. Grüneis, G. H. Booth, G. Kresse, and A. Alavi, “Convergence of many-body wave-function expansions using a plane-wave basis: From homogeneous electron gas to solid state systems,” Phys. Rev. B 86, 035111 (2012a).
  * Booth  _et al._ [2013] G. H. Booth, A. Grüneis, G. Kresse, and A. Alavi, “Towards an exact description of electronic wavefunctions in real solids,” Nature 493, 365–370 (2013).
  * Yang  _et al._ [2014] J. Yang, W. Hu, D. Usvyat, D. Matthews, M. Schütz, and G. K.-L. Chan, “Ab initio determination of the crystalline benzene lattice energy to sub-kilojoule/mole accuracy,” Science 345, 640–643 (2014).
  * McClain  _et al._ [2017] J. McClain, Q. Sun, G. K. L. Chan, and T. C. Berkelbach, “Gaussian-Based Coupled-Cluster Theory for the Ground-State and Band Structure of Solids,” J. Chem. Theory Comput. 13, 1209–1218 (2017).
  * Gruber  _et al._ [2018] T. Gruber, K. Liao, T. Tsatsoulis, F. Hummel, and A. Grüneis, “Applying the Coupled-Cluster Ansatz to Solids and Surfaces in the Thermodynamic Limit,” Phys. Rev. X 8, 021043 (2018).
  * Brandenburg  _et al._ [2019] J. G. Brandenburg, A. Zen, M. Fitzner, B. Ramberger, G. Kresse, T. Tsatsoulis, A. Grüneis, A. Michaelides, and D. Alfè, “Physisorption of Water on Graphene: Subchemical Accuracy from Many-Body Electronic Structure Methods,” J. Phys. Chem. Lett. 10, 358–368 (2019).
  * Tsatsoulis  _et al._ [2017] T. Tsatsoulis, F. Hummel, D. Usvyat, M. Schütz, G. H. Booth, S. S. Binnie, M. J. Gillan, D. Alfè, A. Michaelides, and A. Grüneis, “A comparison between quantum chemistry and quantum Monte Carlo techniques for the adsorption of water on the (001) LiH surface,” J. Chem. Phys. 146, 204108 (2017).
  * Usvyat, Maschio, and Schütz [2018] D. Usvyat, L. Maschio, and M. Schütz, “Periodic and fragment models based on the local correlation approach,” WIRES Comput. Mol. Sci. 8, e1357 (2018).
  * Harl and Kresse [2009] J. Harl and G. Kresse, “Accurate Bulk Properties from Approximate Many-Body Techniques,” Phys. Rev. Lett. 103, 056401 (2009).
  * Grüneis  _et al._ [2009] A. Grüneis, M. Marsman, J. Harl, L. Schimka, and G. Kresse, “Making the random phase approximation to electronic correlation accurate,” J. Chem. Phys. 131 (2009), 10.1063/1.3250347.
  * Shavitt and Bartlett [2009] I. Shavitt and R. J. Bartlett,  _Many-body methods in chemistry and physics: MBPT and coupled-cluster theory_ (Cambridge University Press, Cambridge, 2009).
  * Bartlett [2012] R. J. Bartlett, “Coupled-cluster theory and its equation-of-motion extensions,” WIRES Comput. Mol. Sci. 2, 126–138 (2012).
  * Zhang and Grüneis [2019] I. Y. Zhang and A. Grüneis, “Coupled cluster theory in materials science,” Front. Mater. 6, 123 (2019).
  * Grüneis [2020] A. Grüneis, “Coupled Cluster and Quantum Chemistry Schemes for Solids,” in  _Handbook of Materials Modeling_ (Springer International Publishing, Cham, 2020) pp. 453–468.
  * Stoll, Paulus, and Fulde [2009] H. Stoll, B. Paulus, and P. Fulde, “An incremental coupled-cluster approach to metallic lithium,” Chem. Phys. Lett. 469, 90–93 (2009).
  * Hummel [2018] F. Hummel, “Finite Temperature Coupled Cluster Theories for Extended Systems,” J. Chem. Theory Comput. 14, 6505–6514 (2018).
  * Gell-Mann and Brueckner [1957] M. Gell-Mann and K. A. Brueckner, “Correlation energy of an electron gas at high density,” Phys. Rev. 106, 364–368 (1957).
  * Bohm, Huang, and Pines [1957] D. Bohm, K. Huang, and D. Pines, “Role of Subsidiary Conditions in the Collective Description of Electron Interactions,” Phys. Rev. 107, 71–80 (1957).
  * Freeman [1977] D. L. Freeman, “Coupled-cluster expansion applied to the electron gas: Inclusion of ring and exchange effects,” Phys. Rev. B 15, 5512–5521 (1977).
  * Bartlett and Musiał [2007] R. J. Bartlett and M. Musiał, “Coupled-cluster theory in quantum chemistry,” Rev. Mod. Phys. 79, 291–352 (2007).
  * Shepherd and Grüneis [2013] J. J. Shepherd and A. Grüneis, “Many-body quantum chemistry for the electron gas: Convergent perturbative theories,” Phys. Rev. Lett. 110, 226401 (2013).
  * Liao and Grüneis [2016] K. Liao and A. Grüneis, “Communication: Finite size correction in periodic coupled cluster theory calculations of solids,” J. Chem. Phys. 145, 141102 (2016).
  * Gruber and Grüneis [2018] T. Gruber and A. Grüneis, “Ab initio calculations of carbon and boron nitride allotropes and their structural phase transitions using periodic coupled cluster theory,” Phys. Rev. B 98, 134108 (2018).
  * Dittmer  _et al._ [2019] A. Dittmer, R. Izsák, F. Neese, and D. Maganas, “Accurate Band Gap Predictions of Semiconductors in the Framework of the Similarity Transformed Equation of Motion Coupled Cluster Theory,” Inorg. Chem. 58, 9303–9315 (2019).
  * Gao  _et al._ [2020] Y. Gao, Q. Sun, J. M. Yu, M. Motta, J. McClain, A. F. White, A. J. Minnich, and G. K. L. Chan, “Electronic structure of bulk manganese oxide and nickel oxide from coupled cluster theory,” Phys. Rev. B 101, 165138 (2020).
  * Wang and Berkelbach [2020] X. Wang and T. C. Berkelbach, “Excitons in Solids from Periodic Equation-of-Motion Coupled-Cluster Theory,” J. Chem. Theory Comput. 16, 3095–3103 (2020).
  * Pulkin and Chan [2020] A. Pulkin and G. K. L. Chan, “First-principles coupled cluster theory of the electronic spectrum of transition metal dichalcogenides,” Phys. Rev. B 101, 241113 (2020).
  * Freeman [1978] D. L. Freeman, “Application of the coupled-cluster expansion to the correlation energy of electrons in two-dimensional and quasi-two-dimensional systems,” Solid State Commun. 26, 289–293 (1978).
  * Freeman [1983] D. L. Freeman, “Coupled-cluster summation of the particle-particle ladder diagrams for the two-dimensional electron gas,” J. Phys. C Solid State 16, 711–727 (1983).
  * Shepherd [2016] J. J. Shepherd, “Communication: Convergence of many-body wave-function expansions using a plane-wave basis in the thermodynamic limit,” J. Chem. Phys. 145, 031104 (2016).
  * McClain  _et al._ [2016] J. McClain, J. Lischner, T. Watson, D. A. Matthews, E. Ronca, S. G. Louie, T. C. Berkelbach, and G. K. L. Chan, “Spectral functions of the uniform electron gas via coupled-cluster theory and comparison to the GW and related approximations,” Phys. Rev. B 93, 235139 (2016).
  * Spencer and Thom [2016] J. S. Spencer and A. J. W. Thom, “Developments in stochastic coupled cluster theory: The initiator approximation and application to the uniform electron gas,” J. Chem. Phys. 144, 084108 (2016).
  * Neufeld and Thom [2017] V. A. Neufeld and A. J. Thom, “A study of the dense uniform electron gas with high orders of coupled cluster,” J. Chem. Phys. 147, 194105 (2017).
  * Lewis and Berkelbach [2019] A. M. Lewis and T. C. Berkelbach, “Ab Initio Lifetime and Concomitant Double-Excitation Character of Plasmons at Metallic Densities,” Phys. Rev. Lett. 122, 226402 (2019).
  * White and Kin-Lic Chan [2020] A. F. White and G. Kin-Lic Chan, “Finite-temperature coupled cluster: Efficient implementation and application to prototypical systems,” J. Chem. Phys. 152, 224104 (2020).
  * Shepherd  _et al._ [2012b] J. J. Shepherd, G. Booth, A. Grüneis, and A. Alavi, “Full configuration interaction perspective on the homogeneous electron gas,” Phys. Rev. B 85, 081103 (2012b).
  * Hättig  _et al._ [2012] C. Hättig, W. Klopper, A. Köhn, and D. P. Tew, “Explicitly Correlated Electrons in Molecules,” Chem. Rev. 112, 4–74 (2012).
  * Grüneis  _et al._ [2013] A. Grüneis, J. J. Shepherd, A. Alavi, D. P. Tew, and G. H. Booth, “Explicitly correlated plane waves: Accelerating convergence in periodic wavefunction expansions,” J. Chem. Phys. 139, 084112 (2013).
  * Luo and Alavi [2018] H. Luo and A. Alavi, “Combining the Transcorrelated Method with Full Configuration Interaction Quantum Monte Carlo: Application to the Homogeneous Electron Gas,” J. Chem. Theory Comput. 14, 1403–1411 (2018).
  * [50] K. Liao, T. Schraivogel, H. Luo, D. Kats, and A. Alavi, “Towards efficient and accurate ab initio solutions to periodic systems via transcorrelation and coupled cluster theory,” arXiv:2103.03176 .
  * Fiedler  _et al._ [2017] B. Fiedler, G. Schmitz, C. Hättig, and J. Friedrich, “Combining Accuracy and Efficiency: An Incremental Focal-Point Method Based on Pair Natural Orbitals,” J. Chem. Theory Comput. 13, 6023–6042 (2017).
  * Kumar and Crawford [2017] A. Kumar and T. D. Crawford, “Frozen Virtual Natural Orbitals for Coupled-Cluster Linear-Response Theory,” J. Phys. Chem. A 121, 708–716 (2017).
  * Warden  _et al._ [2020] C. E. Warden, D. G. A. Smith, L. A. Burns, U. Bozkaya, and C. D. Sherrill, “Efficient and automated computation of accurate molecular geometries using focal-point approximations to large-basis coupled-cluster theory,” J. Chem. Phys. 152, 124109 (2020).
  * Irmler  _et al._ [2019] A. Irmler, A. Gallo, F. Hummel, and A. Grüneis, “Duality of Ring and Ladder Diagrams and Its Importance for Many-Electron Perturbation Theories,” Phys. Rev. Lett. 123, 156401 (2019).
  * Irmler and Grüneis [2019] A. Irmler and A. Grüneis, “Particle-particle ladder based basis-set corrections applied to atoms and molecules using coupled-cluster theory,” J. Chem. Phys. 151, 104107 (2019).
  * Lange and Berkelbach [2020] M. F. Lange and T. C. Berkelbach, “Active space approaches combining coupled-cluster and perturbation theory for ground states and excited states,” Mol. Phys. 118, e1808726 (2020).
  * Lau, Knizia, and Berkelbach [2021] B. T. G. Lau, G. Knizia, and T. C. Berkelbach, “Regional Embedding Enables High-Level Quantum Chemistry for Surface Science,” J. Phys. Chem. Lett. 12, 1104–1109 (2021).
  * Kinoshita, Hino, and Bartlett [2005] T. Kinoshita, O. Hino, and R. J. Bartlett, “Coupled-cluster method tailored by configuration interaction,” J. Chem. Phys. 123, 074106 (2005).
  * Piecuch, Kucharski, and Špirko [1999] P. Piecuch, S. A. Kucharski, and V. Špirko, “Coupled-cluster methods with internal and semi-internal triply excited clusters: Vibrational spectrum of the HF molecule,” J. Chem. Phys. 111, 6679–6692 (1999).
  * Kowalski and Piecuch [2001] K. Kowalski and P. Piecuch, “The state-universal multi-reference coupled-cluster theory with perturbative description of core–virtual excitations,” Chem. Phys. Lett. 334, 89–98 (2001).
  * Piecuch [2010] P. Piecuch, “Active-space coupled-cluster methods,” Mol. Phys. 108, 2987–3015 (2010).
  * Shen and Piecuch [2012] J. Shen and P. Piecuch, “Combining active-space coupled-cluster methods with moment energy corrections via the CC(P ;Q) methodology, with benchmark calculations for biradical transition states,” J. Chem. Phys. 136, 144104 (2012).
  * Dutta  _et al._ [2017] A. K. Dutta, M. Nooijen, F. Neese, and R. Izsák, “Automatic active space selection for the similarity transformed equations of motion coupled cluster method,” J. Chem. Phys. 146, 074103 (2017).
  * Nooijen [1999] M. Nooijen, “Combining coupled cluster and perturbation theory,” J. Chem. Phys. 111, 10815–10826 (1999).
  * Bochevarov and Sherrill [2005] A. D. Bochevarov and C. D. Sherrill, “Hybrid correlation models based on active-space partitioning: Correcting second-order Møller-Plesset perturbation theory for bond-breaking reactions,” J. Chem. Phys. 122, 234110 (2005).
  * Myhre, Sánches De Merás, and Koch [2013] R. H. Myhre, A. M. Sánches De Merás, and H. Koch, “The extended CC2 model ECC2,” Mol. Phys. 111, 1109–1118 (2013).
  * Myhre, Sánchez de Merás, and Koch [2014] R. H. Myhre, A. M. J. Sánchez de Merás, and H. Koch, “Multi-level coupled cluster theory,” J. Chem. Phys. 141, 224105 (2014).
  * Aryasetiawan  _et al._ [2004] F. Aryasetiawan, M. Imada, A. Georges, G. Kotliar, S. Biermann, and A. I. Lichtenstein, “Frequency-dependent local interactions and low-energy effective models from electronic structure calculations,” Phys. Rev. B 70, 195104 (2004).
  * Miyake and Aryasetiawan [2008] T. Miyake and F. Aryasetiawan, “Screened Coulomb interaction in the maximally localized Wannier basis,” Phys. Rev. B 77, 085122 (2008).
  * Shepherd, Henderson, and Scuseria [2014] J. J. Shepherd, T. M. Henderson, and G. E. Scuseria, “Coupled cluster channels in the homogeneous electron gas,” J. Chem. Phys. 140, 124102 (2014).
  * Baldereschi [1973] A. Baldereschi, “Mean-value point in the Brillouin zone,” Phys. Rev. B 7, 5212–5215 (1973).
  * Drummond  _et al._ [2008] N. D. Drummond, R. J. Needs, A. Sorouri, and W. M. Foulkes, “Finite-size errors in continuum quantum Monte Carlo calculations,” Phys. Rev. B 78, 125106 (2008).
  * Mihm, McIsaac, and Shepherd [2019] T. N. Mihm, A. R. McIsaac, and J. J. Shepherd, “An optimized twist angle to find the twist-averaged correlation energy applied to the uniform electron gas,” J. Chem. Phys. 150, 191101 (2019).
  * Vosko, Wilk, and Nusair [1980] S. H. Vosko, L. Wilk, and M. Nusair, “Accurate spin-dependent electron liquid correlation energies for local spin density calculations: a critical analysis,” Can. J. Phys. 58, 1200–1211 (1980).
  * Scuseria, Henderson, and Sorensen [2008] G. E. Scuseria, T. M. Henderson, and D. C. Sorensen, “The ground state correlation energy of the random phase approximation from a ring coupled cluster doubles approach,” J. Chem. Phys. 129, 231101 (2008).
  * Schäfer  _et al._ [2021] T. Schäfer, F. Libisch, G. Kresse, and A. Grüneis, “Local embedding of coupled cluster theory into the random phase approximation using plane waves,” J. Chem. Phys. 154, 011101 (2021).
  * Grüneis  _et al._ [2011] A. Grüneis, G. H. Booth, M. Marsman, J. Spencer, A. Alavi, and G. Kresse, “Natural Orbitals for Wave Function Based Correlated Calculations Using a Plane Wave Basis Set,” J. Chem. Theory Comput. 7, 2780–2785 (2011).
  * Sholl [1967] C. A. Sholl, “The calculation of electrostatic energies of metals by plane-wise summation,” Proc. Phys. Soc. 92, 434–445 (1967).