diff --git a/docs/404.html b/docs/404.html index c13cf044..bec63c5d 100644 --- a/docs/404.html +++ b/docs/404.html @@ -1,79 +1,41 @@ - - - - + + + + - Page not found (404) • infer - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - + + + + + + + + + + + + + + - + + + + - - -
-
- +
+ + + + + - - -
+
+
-
- - + - - diff --git a/docs/CONDUCT.html b/docs/CONDUCT.html index 5a380ed8..b3f033d6 100644 --- a/docs/CONDUCT.html +++ b/docs/CONDUCT.html @@ -6,7 +6,7 @@ -Contributor Code of Conduct • infer +Contributor Covenant Code of Conduct • infer @@ -18,10 +18,9 @@ - - - + + @@ -33,7 +32,7 @@ - + @@ -42,17 +41,15 @@ + + + - - - - - - + @@ -66,11 +63,25 @@ + + + + + + + + +
+ +
@@ -153,17 +159,90 @@
-
+
-

As contributors and maintainers of this project, we pledge to respect all people who contribute through reporting issues, posting feature requests, updating documentation, submitting pull requests or patches, and other activities.

-

We are committed to making participation in this project a harassment-free experience for everyone, regardless of level of experience, gender, gender identity and expression, sexual orientation, disability, personal appearance, body size, race, ethnicity, age, or religion.

-

Examples of unacceptable behavior by participants include the use of sexual language or imagery, derogatory comments or personal attacks, trolling, public or private harassment, insults, or other unprofessional conduct.

-

Project maintainers have the right and responsibility to remove, edit, or reject comments, commits, code, wiki edits, issues, and other contributions that are not aligned to this Code of Conduct. Project maintainers who do not follow the Code of Conduct may be removed from the project team.

-

Instances of abusive, harassing, or otherwise unacceptable behavior may be reported by opening an issue or contacting one or more of the project maintainers.

-

This Code of Conduct is adapted from the Contributor Covenant (http://contributor-covenant.org), version 1.0.0, available at http://contributor-covenant.org/version/1/0/0/.

+
+

+Our Pledge

+

We as members, contributors, and leaders pledge to make participation in our community a harassment-free experience for everyone, regardless of age, body size, visible or invisible disability, ethnicity, sex characteristics, gender identity and expression, level of experience, education, socio-economic status, nationality, personal appearance, race, religion, or sexual identity and orientation.

+

We pledge to act and interact in ways that contribute to an open, welcoming, diverse, inclusive, and healthy community.

+
+
+

+Our Standards

+

Examples of behavior that contributes to a positive environment for our community include:

+
    +
  • Demonstrating empathy and kindness toward other people
  • +
  • Being respectful of differing opinions, viewpoints, and experiences
  • +
  • Giving and gracefully accepting constructive feedback
  • +
  • Accepting responsibility and apologizing to those affected by our mistakes, and learning from the experience
  • +
  • Focusing on what is best not just for us as individuals, but for the overall community
  • +
+

Examples of unacceptable behavior include:

+
    +
  • The use of sexualized language or imagery, and sexual attention or advances of any kind
  • +
  • Trolling, insulting or derogatory comments, and personal or political attacks
  • +
  • Public or private harassment
  • +
  • Publishing others’ private information, such as a physical or email address, without their explicit permission
  • +
  • Other conduct which could reasonably be considered inappropriate in a professional setting
  • +
+
+
+

+Enforcement Responsibilities

+

Community leaders are responsible for clarifying and enforcing our standards of acceptable behavior and will take appropriate and fair corrective action in response to any behavior that they deem inappropriate, threatening, offensive, or harmful.

+

Community leaders have the right and responsibility to remove, edit, or reject comments, commits, code, wiki edits, issues, and other contributions that are not aligned to this Code of Conduct, and will communicate reasons for moderation decisions when appropriate.

+
+
+

+Scope

+

This Code of Conduct applies within all community spaces, and also applies when an individual is officially representing the community in public spaces. Examples of representing our community include using an official e-mail address, posting via an official social media account, or acting as an appointed representative at an online or offline event.

+
+
+

+Enforcement

+

Instances of abusive, harassing, or otherwise unacceptable behavior may be reported to the community leaders responsible for enforcement listed as package authors. All complaints will be reviewed and investigated promptly and fairly.

+

All community leaders are obligated to respect the privacy and security of the reporter of any incident.

+
+
+

+Enforcement Guidelines

+

Community leaders will follow these Community Impact Guidelines in determining the consequences for any action they deem in violation of this Code of Conduct:

+
+

+1. Correction

+

Community Impact: Use of inappropriate language or other behavior deemed unprofessional or unwelcome in the community.

+

Consequence: A private, written warning from community leaders, providing clarity around the nature of the violation and an explanation of why the behavior was inappropriate. A public apology may be requested.

+
+
+

+2. Warning

+

Community Impact: A violation through a single incident or series of actions.

+

Consequence: A warning with consequences for continued behavior. No interaction with the people involved, including unsolicited interaction with those enforcing the Code of Conduct, for a specified period of time. This includes avoiding interactions in community spaces as well as external channels like social media. Violating these terms may lead to a temporary or permanent ban.

+
+
+

+3. Temporary Ban

+

Community Impact: A serious violation of community standards, including sustained inappropriate behavior.

+

Consequence: A temporary ban from any sort of interaction or public communication with the community for a specified period of time. No public or private interaction with the people involved, including unsolicited interaction with those enforcing the Code of Conduct, is allowed during this period. Violating these terms may lead to a permanent ban.

+
+
+

+4. Permanent Ban

+

Community Impact: Demonstrating a pattern of violation of community standards, including sustained inappropriate behavior, harassment of an individual, or aggression toward or disparagement of classes of individuals.

+

Consequence: A permanent ban from any sort of public interaction within the community.

+
+
+
+

+Attribution

+

This Code of Conduct is adapted from the Contributor Covenant, version 2.0, available at https://www.contributor-covenant.org/version/2/0/ code_of_conduct.html.

+

Community Impact Guidelines were inspired by Mozilla’s code of conduct enforcement ladder.

+

For answers to common questions about this code of conduct, see the FAQ at https://www.contributor-covenant.org/faq. Translations are available at https:// www.contributor-covenant.org/translations.

+
@@ -179,36 +258,24 @@

Contents

-
- - + diff --git a/docs/CONTRIBUTING.html b/docs/CONTRIBUTING.html index 2111d552..81e5dc92 100644 --- a/docs/CONTRIBUTING.html +++ b/docs/CONTRIBUTING.html @@ -18,10 +18,9 @@ - - - + + @@ -33,7 +32,7 @@ - + @@ -42,13 +41,11 @@ + + + - - - - - @@ -66,11 +63,25 @@ + + + + + + + + +
+ +
@@ -162,23 +168,18 @@

Contributing

Please use the GitHub issues. For any pull request, please link to or open a corresponding issue in GitHub issues. Please ensure that you have notifications turned on and respond to questions, comments or needed changes promptly.

-Tests

+Tests

infer uses testthat for testing. Please try to provide 100% test coverage for any submitted code and always check that existing tests continue to pass. If you are a beginner and need help with writing a test, mention this in the issue and we will try to help.

-

It’s also helpful to run goodpractice::gp() to ensure that lines of code are not over 80 columns and that all lines of code have tests written. Please do so prior to submitting any pull request and fix any suggestions from there. Reach out to us if you need any assistance there too.

-
-
-

-Pull requests

-

Pull requests should be against the develop branch NOT the master branch. You can set this when creating your pull request. Please make a separately named branch to submit. Keep each branch for a complete specific issue.

+

It’s also helpful to run goodpractice::gp() to ensure that lines of code are not over 80 characters and that all lines of code have tests written. Please do so prior to submitting any pull request and fix any suggestions from there. Reach out to us if you need any assistance there too.

-Code style

-

Please use snake case (such as rep_sample_n) for function names. Besides that, in general follow the tidyverse style for R.

+Code style +

Please use snake case (such as rep_sample_n) for function names. Besides that, in general follow the tidyverse style for R.

-Code of Conduct

+Code of Conduct

When contributing to the infer package you must follow the code of conduct defined in CONDUCT.

@@ -196,36 +197,24 @@

Contents

- - + diff --git a/docs/articles/anova.html b/docs/articles/anova.html index 04f61984..a51fa2fa 100644 --- a/docs/articles/anova.html +++ b/docs/articles/anova.html @@ -12,21 +12,30 @@ - + - - - - + + + + + + + +
@@ -105,13 +109,13 @@ -
+
@@ -120,73 +124,87 @@

Tidy ANOVA (Analysis of Variance) with infer

In this vignette, we’ll walk through conducting an analysis of variance (ANOVA) test using infer. ANOVAs are used to analyze differences in group means.

Throughout this vignette, we’ll make use of the gss dataset supplied by infer, which contains a sample of data from the General Social Survey. See ?gss for more information on the variables included and their source. Note that this data (and our examples on it) are for demonstration purposes only, and will not necessarily provide accurate estimates unless weighted properly. For these examples, let’s suppose that this dataset is a representative sample of a population we want to learn about: American adults. The data looks like this:

-
dplyr::glimpse(gss)
+
+dplyr::glimpse(gss)
## Rows: 500
 ## Columns: 11
-## $ year    <dbl> 2014, 1994, 1998, 1996, 1994, 1996, 1990, 2016, 2000, 1998, 2…
-## $ age     <dbl> 36, 34, 24, 42, 31, 32, 48, 36, 30, 33, 21, 30, 38, 49, 25, 5…
-## $ sex     <fct> male, female, male, male, male, female, female, female, femal…
-## $ college <fct> degree, no degree, degree, no degree, degree, no degree, no d…
-## $ partyid <fct> ind, rep, ind, ind, rep, rep, dem, ind, rep, dem, dem, ind, d…
-## $ hompop  <dbl> 3, 4, 1, 4, 2, 4, 2, 1, 5, 2, 4, 3, 4, 4, 2, 2, 3, 2, 1, 2, 5…
-## $ hours   <dbl> 50, 31, 40, 40, 40, 53, 32, 20, 40, 40, 23, 52, 38, 72, 48, 4…
-## $ income  <ord> $25000 or more, $20000 - 24999, $25000 or more, $25000 or mor…
-## $ class   <fct> middle class, working class, working class, working class, mi…
-## $ finrela <fct> below average, below average, below average, above average, a…
-## $ weight  <dbl> 0.8960, 1.0825, 0.5501, 1.0864, 1.0825, 1.0864, 1.0627, 0.478…
+## $ year <dbl> 2014, 1994, 1998, 1996, 1994, 1996, 1990, 2016, 2000, 1998, 20… +## $ age <dbl> 36, 34, 24, 42, 31, 32, 48, 36, 30, 33, 21, 30, 38, 49, 25, 56… +## $ sex <fct> male, female, male, male, male, female, female, female, female… +## $ college <fct> degree, no degree, degree, no degree, degree, no degree, no de… +## $ partyid <fct> ind, rep, ind, ind, rep, rep, dem, ind, rep, dem, dem, ind, de… +## $ hompop <dbl> 3, 4, 1, 4, 2, 4, 2, 1, 5, 2, 4, 3, 4, 4, 2, 2, 3, 2, 1, 2, 5,… +## $ hours <dbl> 50, 31, 40, 40, 40, 53, 32, 20, 40, 40, 23, 52, 38, 72, 48, 40… +## $ income <ord> $25000 or more, $20000 - 24999, $25000 or more, $25000 or more… +## $ class <fct> middle class, working class, working class, working class, mid… +## $ finrela <fct> below average, below average, below average, above average, ab… +## $ weight <dbl> 0.8960, 1.0825, 0.5501, 1.0864, 1.0825, 1.0864, 1.0627, 0.4785…

To carry out an ANOVA, we’ll examine the association between age and political party affiliation in the United States. The age variable is a numerical variable measuring the respondents’ age at the time that the survey was taken, and partyid is a factor variable with unique values ind, rep, dem, other.

This is what the relationship looks like in the observed data:

If there were no relationship, we would expect to see the each of these boxplots lining up along the y-axis. It looks like the average age of democrats and republicans seems to be a bit larger than independent and other American voters. Is this difference just random noise, though?

First, to calculate the observed statistic, we can use specify() and calculate().

-
# calculate the observed statistic
-observed_f_statistic <- gss %>%
-  specify(age ~ partyid) %>%
-  calculate(stat = "F")
+
+# calculate the observed statistic
+observed_f_statistic <- gss %>%
+  specify(age ~ partyid) %>%
+  hypothesize(null = "independence") %>%
+  calculate(stat = "F")

The observed \(F\) statistic is 2.4842. Now, we want to compare this statistic to a null distribution, generated under the assumption that age and political party affiliation are not actually related, to get a sense of how likely it would be for us to see this observed statistic if there were actually no association between the two variables.

-

We can generate the null distribution using randomization. The randomization approach permutes the response and explanatory variables, so that each person’s educational attainment is matched up with a random income from the sample in order to break up any association between the two.

-
# generate the null distribution using randomization
-null_distribution <- gss %>%
-  specify(age ~ partyid) %>%
-  hypothesize(null = "independence") %>%
-  generate(reps = 1000, type = "permute") %>%
-  calculate(stat = "F")
-

Note that, in the line specify(age ~ partyid) above, we could use the equivalent syntax specify(response = age, explanatory = partyid).

+

We can generate an approximation of the null distribution using randomization. The randomization approach permutes the response and explanatory variables, so that each person’s party affiliation is matched up with a random age from the sample in order to break up any association between the two.

+
+# generate the null distribution using randomization
+null_dist <- gss %>%
+  specify(age ~ partyid) %>%
+  hypothesize(null = "independence") %>%
+  generate(reps = 1000, type = "permute") %>%
+  calculate(stat = "F")
+

Note that, in the line specify(age ~ partyid) above, we could use the equivalent syntax specify(response = age, explanatory = partyid).

To get a sense for what this distribution looks like, and where our observed statistic falls, we can use visualize():

-
# visualize the null distribution and test statistic!
-null_distribution %>%
-  visualize() +
-  shade_p_value(observed_f_statistic,
-                direction = "greater")
+
+# visualize the null distribution and test statistic!
+null_dist %>%
+  visualize() + 
+  shade_p_value(observed_f_statistic,
+                direction = "greater")

-

We could also visualize the observed statistic against the theoretical null distribution. Note that we skip the generate() and calculate() steps when using the theoretical approach, and that we now need to provide method = "theoretical" to visualize().

-
# visualize the theoretical null distribution and test statistic!
-gss %>%
-  specify(age ~ partyid) %>%
-  hypothesize(null = "independence") %>%
-  visualize(method = "theoretical") +
-  shade_p_value(observed_f_statistic,
-                direction = "greater")
+

We could also visualize the observed statistic against the theoretical null distribution. To do so, use the assume() verb to define a theoretical null distribution and then pass it to visualize() like a null distribution outputted from generate() and calculate().

+
+# visualize the theoretical null distribution and test statistic!
+null_dist_theory <- gss %>%
+  specify(age ~ partyid) %>%
+  assume(distribution = "F")
+
+visualize(null_dist_theory) +
+  shade_p_value(observed_f_statistic,
+                direction = "greater")

To visualize both the randomization-based and theoretical null distributions to get a sense of how the two relate, we can pipe the randomization-based null distribution into visualize(), and then further provide method = "both" to visualize().

-
# visualize both null distributions and the test statistic!
-null_distribution %>%
-  visualize(method = "both") +
-  shade_p_value(observed_f_statistic,
-                direction = "greater")
+
+# visualize both null distributions and the test statistic!
+null_dist %>%
+  visualize(method = "both") + 
+  shade_p_value(observed_f_statistic,
+                direction = "greater")

-

Either way, it looks like our observed test statistic would be really unlikely if there were actually no association between age and political party affiliation. More exactly, we can calculate the p-value:

-
# calculate the p value from the observed statistic and null distribution
-p_value <- null_distribution %>%
-  get_p_value(obs_stat = observed_f_statistic,
-              direction = "greater")
-
-p_value
-
## # A tibble: 1 x 1
+

Either way, it looks like our observed test statistic would be quite unlikely if there were actually no association between age and political party affiliation. More exactly, we can approximate the p-value from the randomization-based approximation to the null distribution:

+
+# calculate the p value from the observed statistic and null distribution
+p_value <- null_dist %>%
+  get_p_value(obs_stat = observed_f_statistic,
+              direction = "greater")
+
+p_value
+
## # A tibble: 1 × 1
 ##   p_value
 ##     <dbl>
-## 1   0.056
-

Thus, if there were really no relationship between age and political party affiliation, the probability that we would see a statistic as or more extreme than 2.4842 is approximately 0.056.

+## 1 0.055
+

Thus, if there were really no relationship between age and political party affiliation, our approximation of the probability that we would see a statistic as or more extreme than 2.4842 is approximately 0.055.

+

To calculate the p-value using the true \(F\) distribution, we can use the pf function from base R. This function allows us to situate the test statistic we calculated previously in the \(F\) distribution with the appropriate degrees of freedom.

+
+pf(observed_f_statistic$stat, 3, 496, lower.tail = FALSE)
+
## [1] 0.06006
+

Note that, while the observed statistic stays the same, the resulting p-value differs slightly between these two approaches since the randomization-based empirical \(F\) distribution is an approximation of the true \(F\) distribution.

The package currently does not supply a wrapper for tidy ANOVA tests.

@@ -198,32 +216,24 @@

Tidy ANOVA (Analysis of Variance) with infer

-
- - + + + + + diff --git a/docs/articles/anova_files/figure-html/plot-f-1.png b/docs/articles/anova_files/figure-html/plot-f-1.png index f0b81d3d..b96d15af 100644 Binary files a/docs/articles/anova_files/figure-html/plot-f-1.png and b/docs/articles/anova_files/figure-html/plot-f-1.png differ diff --git a/docs/articles/anova_files/figure-html/visualize-f-1.png b/docs/articles/anova_files/figure-html/visualize-f-1.png index 0149145b..51f82560 100644 Binary files a/docs/articles/anova_files/figure-html/visualize-f-1.png and b/docs/articles/anova_files/figure-html/visualize-f-1.png differ diff --git a/docs/articles/anova_files/figure-html/visualize-f-theor-1.png b/docs/articles/anova_files/figure-html/visualize-f-theor-1.png index 0059ae73..e89c374c 100644 Binary files a/docs/articles/anova_files/figure-html/visualize-f-theor-1.png and b/docs/articles/anova_files/figure-html/visualize-f-theor-1.png differ diff --git a/docs/articles/anova_files/figure-html/visualize-indep-both-1.png b/docs/articles/anova_files/figure-html/visualize-indep-both-1.png index a77ac2a5..10cfbe75 100644 Binary files a/docs/articles/anova_files/figure-html/visualize-indep-both-1.png and b/docs/articles/anova_files/figure-html/visualize-indep-both-1.png differ diff --git a/docs/articles/anova_files/header-attrs-2.8/header-attrs.js b/docs/articles/anova_files/header-attrs-2.8/header-attrs.js new file mode 100644 index 00000000..dd57d92e --- /dev/null +++ b/docs/articles/anova_files/header-attrs-2.8/header-attrs.js @@ -0,0 +1,12 @@ +// Pandoc 2.9 adds attributes on both header and div. We remove the former (to +// be compatible with the behavior of Pandoc < 2.8). +document.addEventListener('DOMContentLoaded', function(e) { + var hs = document.querySelectorAll("div.section[class*='level'] > :first-child"); + var i, h, a; + for (i = 0; i < hs.length; i++) { + h = hs[i]; + if (!/^h[1-6]$/i.test(h.tagName)) continue; // it should be a header h1-h6 + a = h.attributes; + while (a.length > 0) h.removeAttribute(a[0].name); + } +}); diff --git a/docs/articles/chi_squared.html b/docs/articles/chi_squared.html index 38061d1f..5b18877e 100644 --- a/docs/articles/chi_squared.html +++ b/docs/articles/chi_squared.html @@ -12,21 +12,30 @@ - + - - - - + + + + + + + +
@@ -105,13 +109,13 @@ -
+
@@ -120,153 +124,175 @@

Tidy Chi-Squared Tests with infer

-Introduction

+Introduction

In this vignette, we’ll walk through conducting a \(\chi^2\) (chi-squared) test of independence and a chi-squared goodness of fit test using infer. We’ll start out with a chi-squared test of independence, which can be used to test the association between two categorical variables. Then, we’ll move on to a chi-squared goodness of fit test, which tests how well the distribution of one categorical variable can be approximated by some theoretical distribution.

Throughout this vignette, we’ll make use of the gss dataset supplied by infer, which contains a sample of data from the General Social Survey. See ?gss for more information on the variables included and their source. Note that this data (and our examples on it) are for demonstration purposes only, and will not necessarily provide accurate estimates unless weighted properly. For these examples, let’s suppose that this dataset is a representative sample of a population we want to learn about: American adults. The data looks like this:

-
dplyr::glimpse(gss)
+
+dplyr::glimpse(gss)
## Rows: 500
 ## Columns: 11
-## $ year    <dbl> 2014, 1994, 1998, 1996, 1994, 1996, 1990, 2016, 2000, 1998, 2…
-## $ age     <dbl> 36, 34, 24, 42, 31, 32, 48, 36, 30, 33, 21, 30, 38, 49, 25, 5…
-## $ sex     <fct> male, female, male, male, male, female, female, female, femal…
-## $ college <fct> degree, no degree, degree, no degree, degree, no degree, no d…
-## $ partyid <fct> ind, rep, ind, ind, rep, rep, dem, ind, rep, dem, dem, ind, d…
-## $ hompop  <dbl> 3, 4, 1, 4, 2, 4, 2, 1, 5, 2, 4, 3, 4, 4, 2, 2, 3, 2, 1, 2, 5…
-## $ hours   <dbl> 50, 31, 40, 40, 40, 53, 32, 20, 40, 40, 23, 52, 38, 72, 48, 4…
-## $ income  <ord> $25000 or more, $20000 - 24999, $25000 or more, $25000 or mor…
-## $ class   <fct> middle class, working class, working class, working class, mi…
-## $ finrela <fct> below average, below average, below average, above average, a…
-## $ weight  <dbl> 0.8960, 1.0825, 0.5501, 1.0864, 1.0825, 1.0864, 1.0627, 0.478…
+## $ year <dbl> 2014, 1994, 1998, 1996, 1994, 1996, 1990, 2016, 2000, 1998, 20… +## $ age <dbl> 36, 34, 24, 42, 31, 32, 48, 36, 30, 33, 21, 30, 38, 49, 25, 56… +## $ sex <fct> male, female, male, male, male, female, female, female, female… +## $ college <fct> degree, no degree, degree, no degree, degree, no degree, no de… +## $ partyid <fct> ind, rep, ind, ind, rep, rep, dem, ind, rep, dem, dem, ind, de… +## $ hompop <dbl> 3, 4, 1, 4, 2, 4, 2, 1, 5, 2, 4, 3, 4, 4, 2, 2, 3, 2, 1, 2, 5,… +## $ hours <dbl> 50, 31, 40, 40, 40, 53, 32, 20, 40, 40, 23, 52, 38, 72, 48, 40… +## $ income <ord> $25000 or more, $20000 - 24999, $25000 or more, $25000 or more… +## $ class <fct> middle class, working class, working class, working class, mid… +## $ finrela <fct> below average, below average, below average, above average, ab… +## $ weight <dbl> 0.8960, 1.0825, 0.5501, 1.0864, 1.0825, 1.0864, 1.0627, 0.4785…

-Test of Independence

+Test of Independence

To carry out a chi-squared test of independence, we’ll examine the association between income and educational attainment in the United States. college is a categorical variable with values degree and no degree, indicating whether or not the respondent has a college degree (including community college), and finrela gives the respondent’s self-identification of family income—either far below average, below average, average, above average, far above average, or DK (don’t know).

This is what the relationship looks like in the sample data:

If there were no relationship, we would expect to see the purple bars reaching to the same height, regardless of income class. Are the differences we see here, though, just due to random noise?

First, to calculate the observed statistic, we can use specify() and calculate().

-
# calculate the observed statistic
-observed_indep_statistic <- gss %>%
-  specify(college ~ finrela) %>%
-  calculate(stat = "Chisq")
+
+# calculate the observed statistic
+observed_indep_statistic <- gss %>%
+  specify(college ~ finrela) %>%
+  hypothesize(null = "independence") %>%
+  calculate(stat = "Chisq")

The observed \(\chi^2\) statistic is 30.6825. Now, we want to compare this statistic to a null distribution, generated under the assumption that these variables are not actually related, to get a sense of how likely it would be for us to see this observed statistic if there were actually no association between education and income.

-

We can generate the null distribution in one of two ways—using randomization or theory-based methods. The randomization approach permutes the response and explanatory variables, so that each person’s educational attainment is matched up with a random income from the sample in order to break up any association between the two.

-
# generate the null distribution using randomization
-null_distribution_simulated <- gss %>%
-  specify(college ~ finrela) %>%
-  hypothesize(null = "independence") %>%
-  generate(reps = 1000, type = "permute") %>%
-  calculate(stat = "Chisq")
-

Note that, in the line specify(college ~ finrela) above, we could use the equivalent syntax specify(response = college, explanatory = finrela). The same goes in the code below, which generates the null distribution using theory-based methods instead of randomization.

-
# generate the null distribution by theoretical approximation
-null_distribution_theoretical <- gss %>%
-  specify(college ~ finrela) %>%
-  hypothesize(null = "independence") %>%
-  # note that we skip the generation step here!
-  calculate(stat = "Chisq")
+

We can generate the null distribution in one of two ways—using randomization or theory-based methods. The randomization approach approximates the null distribution by permuting the response and explanatory variables, so that each person’s educational attainment is matched up with a random income from the sample in order to break up any association between the two.

+
+# generate the null distribution using randomization
+null_dist_sim <- gss %>%
+  specify(college ~ finrela) %>%
+  hypothesize(null = "independence") %>%
+  generate(reps = 1000, type = "permute") %>%
+  calculate(stat = "Chisq")
+

Note that, in the line specify(college ~ finrela) above, we could use the equivalent syntax specify(response = college, explanatory = finrela). The same goes in the code below, which generates the null distribution using theory-based methods instead of randomization.

+
+# generate the null distribution by theoretical approximation
+null_dist_theory <- gss %>%
+  specify(college ~ finrela) %>%
+  assume(distribution = "Chisq")

To get a sense for what these distributions look like, and where our observed statistic falls, we can use visualize():

-
# visualize the null distribution and test statistic!
-null_distribution_simulated %>%
-  visualize() +
-  shade_p_value(observed_indep_statistic,
-                direction = "greater")
+
+# visualize the null distribution and test statistic!
+null_dist_sim %>%
+  visualize() + 
+  shade_p_value(observed_indep_statistic,
+                direction = "greater")

-

We could also visualize the observed statistic against the theoretical null distribution. Note that we skip the generate() and calculate() steps when using the theoretical approach, and that we now need to provide method = "theoretical" to visualize().

-
# visualize the theoretical null distribution and test statistic!
-gss %>%
-  specify(college ~ finrela) %>%
-  hypothesize(null = "independence") %>%
-  visualize(method = "theoretical") +
-  shade_p_value(observed_indep_statistic,
-                direction = "greater")
+

We could also visualize the observed statistic against the theoretical null distribution. To do so, use the assume() verb to define a theoretical null distribution and then pass it to visualize() like a null distribution outputted from generate() and calculate().

+
+# visualize the theoretical null distribution and test statistic!
+gss %>%
+  specify(college ~ finrela) %>%
+  assume(distribution = "Chisq") %>%
+  visualize() + 
+  shade_p_value(observed_indep_statistic,
+                direction = "greater")

To visualize both the randomization-based and theoretical null distributions to get a sense of how the two relate, we can pipe the randomization-based null distribution into visualize(), and further provide method = "both".

-
# visualize both null distributions and the test statistic!
-null_distribution_simulated %>%
-  visualize(method = "both") +
-  shade_p_value(observed_indep_statistic,
-                direction = "greater")
+
+# visualize both null distributions and the test statistic!
+null_dist_sim %>%
+  visualize(method = "both") + 
+  shade_p_value(observed_indep_statistic,
+                direction = "greater")

-

Either way, it looks like our observed test statistic would be really unlikely if there were actually no association between education and income. More exactly, we can calculate the p-value:

-
# calculate the p value from the observed statistic and null distribution
-p_value_independence <- null_distribution_simulated %>%
-  get_p_value(obs_stat = observed_indep_statistic,
-              direction = "greater")
+

Either way, it looks like our observed test statistic would be quite unlikely if there were actually no association between education and income. More exactly, we can approximate the p-value with get_p_value:

+
+# calculate the p value from the observed statistic and null distribution
+p_value_independence <- null_dist_sim %>%
+  get_p_value(obs_stat = observed_indep_statistic,
+              direction = "greater")
 
-p_value_independence
-
## # A tibble: 1 x 1
+p_value_independence
+
## # A tibble: 1 × 1
 ##   p_value
 ##     <dbl>
 ## 1       0
-

Thus, if there were really no relationship between education and income, the probability that we would see a statistic as or more extreme than 30.6825 is approximately 0.

-

Note that, equivalently to the steps shown above, the package supplies a wrapper function, chisq_test, to carry out Chi-Squared tests of independence on tidy data. The syntax goes like this:

-
chisq_test(gss, college ~ finrela)
-
## # A tibble: 1 x 3
+

Thus, if there were really no relationship between education and income, our approximation of the probability that we would see a statistic as or more extreme than 30.6825 is approximately 0.

+

To calculate the p-value using the true \(\chi^2\) distribution, we can use the pchisq function from base R. This function allows us to situate the test statistic we calculated previously in the \(\chi^2\) distribution with the appropriate degrees of freedom.

+
+pchisq(observed_indep_statistic$stat, 5, lower.tail = FALSE)
+
## X-squared 
+## 1.082e-05
+

Note that, equivalently to the theory-based approach shown above, the package supplies a wrapper function, chisq_test, to carry out Chi-Squared tests of independence on tidy data. The syntax goes like this:

+
+chisq_test(gss, college ~ finrela)
+
## # A tibble: 1 × 3
 ##   statistic chisq_df   p_value
 ##       <dbl>    <int>     <dbl>
 ## 1      30.7        5 0.0000108

-Goodness of Fit

+Goodness of Fit

Now, moving on to a chi-squared goodness of fit test, we’ll take a look at the self-identified income class of our survey respondents. Suppose our null hypothesis is that finrela follows a uniform distribution (i.e. there’s actually an equal number of people that describe their income as far below average, below average, average, above average, far above average, or that don’t know their income.) The graph below represents this hypothesis:

It seems like a uniform distribution may not be the most appropriate description of the data–many more people describe their income as average than than any of the other options. Lets now test whether this difference in distributions is statistically significant.

First, to carry out this hypothesis test, we would calculate our observed statistic.

-
# calculating the null distribution
-observed_gof_statistic <- gss %>%
-  specify(response = finrela) %>%
-  hypothesize(null = "point",
-              p = c("far below average" = 1/6,
-                    "below average" = 1/6,
-                    "average" = 1/6,
-                    "above average" = 1/6,
-                    "far above average" = 1/6,
-                    "DK" = 1/6)) %>%
-  calculate(stat = "Chisq")
+
+# calculating the null distribution
+observed_gof_statistic <- gss %>%
+  specify(response = finrela) %>%
+  hypothesize(null = "point",
+              p = c("far below average" = 1/6,
+                    "below average" = 1/6,
+                    "average" = 1/6,
+                    "above average" = 1/6,
+                    "far above average" = 1/6,
+                    "DK" = 1/6)) %>%
+  calculate(stat = "Chisq")

The observed statistic is 487.984. Now, generating a null distribution, by just dropping in a call to generate():

-
# generating a null distribution, assuming each income class is equally likely
-null_distribution_gof <- gss %>%
-  specify(response = finrela) %>%
-  hypothesize(null = "point",
-              p = c("far below average" = 1/6,
-                    "below average" = 1/6,
-                    "average" = 1/6,
-                    "above average" = 1/6,
-                    "far above average" = 1/6,
-                    "DK" = 1/6)) %>%
-  generate(reps = 1000, type = "simulate") %>%
-  calculate(stat = "Chisq")
+
+# generating a null distribution, assuming each income class is equally likely
+null_dist_gof <- gss %>%
+  specify(response = finrela) %>%
+  hypothesize(null = "point",
+              p = c("far below average" = 1/6,
+                    "below average" = 1/6,
+                    "average" = 1/6,
+                    "above average" = 1/6,
+                    "far above average" = 1/6,
+                    "DK" = 1/6)) %>%
+  generate(reps = 1000, type = "draw") %>%
+  calculate(stat = "Chisq")

Again, to get a sense for what these distributions look like, and where our observed statistic falls, we can use visualize():

-
# visualize the null distribution and test statistic!
-null_distribution_gof %>%
-  visualize() +
-  shade_p_value(observed_gof_statistic,
-                direction = "greater")
+
+# visualize the null distribution and test statistic!
+null_dist_gof %>%
+  visualize() + 
+  shade_p_value(observed_gof_statistic,
+                direction = "greater")

-

This statistic seems like it would be really unlikely if income class self-identification actually followed a uniform distribution! How unlikely, though? Calculating the p-value:

-
# calculate the p-value
-p_value_gof <- null_distribution_gof %>%
-  get_p_value(observed_gof_statistic,
-              direction = "greater")
+

This statistic seems like it would be quite unlikely if income class self-identification actually followed a uniform distribution! How unlikely, though? Calculating the p-value:

+
+# calculate the p-value
+p_value_gof <- null_dist_gof %>%
+  get_p_value(observed_gof_statistic,
+              direction = "greater")
 
-p_value_gof
-
## # A tibble: 1 x 1
+p_value_gof
+
## # A tibble: 1 × 1
 ##   p_value
 ##     <dbl>
 ## 1       0
-

Thus, if each self-identified income class was equally likely to occur, the probability that we would see a distribution like the one we did is approximately 0.

-

Again, equivalently to the steps shown above, the package supplies a wrapper function, chisq_test, to carry out Chi-Squared goodness of fit tests on tidy data. The syntax goes like this:

-
chisq_test(gss,
-           response = finrela,
-           p = c("far below average" = 1/6,
-                    "below average" = 1/6,
-                    "average" = 1/6,
-                    "above average" = 1/6,
-                    "far above average" = 1/6,
-                    "DK" = 1/6))
-
## # A tibble: 1 x 3
+

Thus, if each self-identified income class was equally likely to occur, our approximation of the probability that we would see a distribution like the one we did is approximately 0.

+

To calculate the p-value using the true \(\chi^2\) distribution, we can use the pchisq function from base R. This function allows us to situate the test statistic we calculated previously in the \(\chi^2\) distribution with the appropriate degrees of freedom.

+
+pchisq(observed_gof_statistic$stat, 5, lower.tail = FALSE)
+
## [1] 3.131e-103
+

Again, equivalently to the theory-based approach shown above, the package supplies a wrapper function, chisq_test, to carry out Chi-Squared goodness of fit tests on tidy data. The syntax goes like this:

+
+chisq_test(gss, 
+           response = finrela,
+           p = c("far below average" = 1/6,
+                    "below average" = 1/6,
+                    "average" = 1/6,
+                    "above average" = 1/6,
+                    "far above average" = 1/6,
+                    "DK" = 1/6))
+
## # A tibble: 1 × 3
 ##   statistic chisq_df   p_value
 ##       <dbl>    <dbl>     <dbl>
 ## 1      488.        5 3.13e-103
@@ -281,32 +307,24 @@

-
- - + + + + + diff --git a/docs/articles/chi_squared_files/figure-html/gof-plot-1.png b/docs/articles/chi_squared_files/figure-html/gof-plot-1.png index 3799073b..9417e615 100644 Binary files a/docs/articles/chi_squared_files/figure-html/gof-plot-1.png and b/docs/articles/chi_squared_files/figure-html/gof-plot-1.png differ diff --git a/docs/articles/chi_squared_files/figure-html/plot-indep-1.png b/docs/articles/chi_squared_files/figure-html/plot-indep-1.png index 4f216984..0d01d1e6 100644 Binary files a/docs/articles/chi_squared_files/figure-html/plot-indep-1.png and b/docs/articles/chi_squared_files/figure-html/plot-indep-1.png differ diff --git a/docs/articles/chi_squared_files/figure-html/visualize-indep-1.png b/docs/articles/chi_squared_files/figure-html/visualize-indep-1.png index 8c556bd5..6c58690c 100644 Binary files a/docs/articles/chi_squared_files/figure-html/visualize-indep-1.png and b/docs/articles/chi_squared_files/figure-html/visualize-indep-1.png differ diff --git a/docs/articles/chi_squared_files/figure-html/visualize-indep-both-1.png b/docs/articles/chi_squared_files/figure-html/visualize-indep-both-1.png index 34984a95..ae2ac078 100644 Binary files a/docs/articles/chi_squared_files/figure-html/visualize-indep-both-1.png and b/docs/articles/chi_squared_files/figure-html/visualize-indep-both-1.png differ diff --git a/docs/articles/chi_squared_files/figure-html/visualize-indep-gof-1.png b/docs/articles/chi_squared_files/figure-html/visualize-indep-gof-1.png index e6948ede..c2a635bb 100644 Binary files a/docs/articles/chi_squared_files/figure-html/visualize-indep-gof-1.png and b/docs/articles/chi_squared_files/figure-html/visualize-indep-gof-1.png differ diff --git a/docs/articles/chi_squared_files/figure-html/visualize-indep-theor-1.png b/docs/articles/chi_squared_files/figure-html/visualize-indep-theor-1.png index 1aeca520..7105fb58 100644 Binary files a/docs/articles/chi_squared_files/figure-html/visualize-indep-theor-1.png and b/docs/articles/chi_squared_files/figure-html/visualize-indep-theor-1.png differ diff --git a/docs/articles/chi_squared_files/header-attrs-2.8/header-attrs.js b/docs/articles/chi_squared_files/header-attrs-2.8/header-attrs.js new file mode 100644 index 00000000..dd57d92e --- /dev/null +++ b/docs/articles/chi_squared_files/header-attrs-2.8/header-attrs.js @@ -0,0 +1,12 @@ +// Pandoc 2.9 adds attributes on both header and div. We remove the former (to +// be compatible with the behavior of Pandoc < 2.8). +document.addEventListener('DOMContentLoaded', function(e) { + var hs = document.querySelectorAll("div.section[class*='level'] > :first-child"); + var i, h, a; + for (i = 0; i < hs.length; i++) { + h = hs[i]; + if (!/^h[1-6]$/i.test(h.tagName)) continue; // it should be a header h1-h6 + a = h.attributes; + while (a.length > 0) h.removeAttribute(a[0].name); + } +}); diff --git a/docs/articles/index.html b/docs/articles/index.html index 460558b7..5909003b 100644 --- a/docs/articles/index.html +++ b/docs/articles/index.html @@ -18,10 +18,9 @@ - - - + + @@ -33,7 +32,7 @@ - + @@ -42,13 +41,11 @@ + + + - - - - - @@ -66,11 +63,25 @@ + + + + + + + +
+ +

@@ -178,36 +184,24 @@

All vignettes

-
- - + diff --git a/docs/articles/infer.html b/docs/articles/infer.html index 6a18495f..c0f9074d 100644 --- a/docs/articles/infer.html +++ b/docs/articles/infer.html @@ -12,21 +12,30 @@ - + - - - - + + + + + + + +
@@ -105,13 +109,13 @@ -
+
@@ -120,7 +124,7 @@

Getting to Know infer

-Introduction

+Introduction

infer implements an expressive grammar to perform statistical inference that coheres with the tidyverse design framework. Rather than providing methods for specific statistical tests, this package consolidates the principles that are shared among common hypothesis tests into a set of 4 main verbs (functions), supplemented with many utilities to visualize and extract value from their outputs.

Regardless of which hypothesis test we’re using, we’re still asking the same kind of question: is the effect/difference in our observed data real, or due to chance? To answer this question, we start by assuming that the observed data came from some world where “nothing is going on” (i.e. the observed effect was simply due to random chance), and call this assumption our null hypothesis. (In reality, we might not believe in the null hypothesis at all—the null hypothesis is in opposition to the alternate hypothesis, which supposes that the effect present in the observed data is actually due to the fact that “something is going on.”) We then calculate a test statistic from our data that describes the observed effect. We can use this test statistic to calculate a p-value, giving the probability that our observed data could come about if the null hypothesis was true. If this probability is below some pre-defined significance level \(\alpha\), then we can reject our null hypothesis.

The workflow of this package is designed around this idea. Starting out with some dataset,

@@ -135,34 +139,36 @@

calculate() allows you to calculate a distribution of statistics from the generated data to form the null distribution.

Throughout this vignette, we make use of gss, a dataset supplied by infer containing a sample of 500 observations of 11 variables from the General Social Survey.

-
# load in the dataset
-data(gss)
+
+# load in the dataset
+data(gss)
 
 # take a look at its structure
-dplyr::glimpse(gss)
+dplyr::glimpse(gss)
## Rows: 500
 ## Columns: 11
-## $ year    <dbl> 2014, 1994, 1998, 1996, 1994, 1996, 1990, 2016, 2000, 1998, 2…
-## $ age     <dbl> 36, 34, 24, 42, 31, 32, 48, 36, 30, 33, 21, 30, 38, 49, 25, 5…
-## $ sex     <fct> male, female, male, male, male, female, female, female, femal…
-## $ college <fct> degree, no degree, degree, no degree, degree, no degree, no d…
-## $ partyid <fct> ind, rep, ind, ind, rep, rep, dem, ind, rep, dem, dem, ind, d…
-## $ hompop  <dbl> 3, 4, 1, 4, 2, 4, 2, 1, 5, 2, 4, 3, 4, 4, 2, 2, 3, 2, 1, 2, 5…
-## $ hours   <dbl> 50, 31, 40, 40, 40, 53, 32, 20, 40, 40, 23, 52, 38, 72, 48, 4…
-## $ income  <ord> $25000 or more, $20000 - 24999, $25000 or more, $25000 or mor…
-## $ class   <fct> middle class, working class, working class, working class, mi…
-## $ finrela <fct> below average, below average, below average, above average, a…
-## $ weight  <dbl> 0.8960, 1.0825, 0.5501, 1.0864, 1.0825, 1.0864, 1.0627, 0.478…
+## $ year <dbl> 2014, 1994, 1998, 1996, 1994, 1996, 1990, 2016, 2000, 1998, 20… +## $ age <dbl> 36, 34, 24, 42, 31, 32, 48, 36, 30, 33, 21, 30, 38, 49, 25, 56… +## $ sex <fct> male, female, male, male, male, female, female, female, female… +## $ college <fct> degree, no degree, degree, no degree, degree, no degree, no de… +## $ partyid <fct> ind, rep, ind, ind, rep, rep, dem, ind, rep, dem, dem, ind, de… +## $ hompop <dbl> 3, 4, 1, 4, 2, 4, 2, 1, 5, 2, 4, 3, 4, 4, 2, 2, 3, 2, 1, 2, 5,… +## $ hours <dbl> 50, 31, 40, 40, 40, 53, 32, 20, 40, 40, 23, 52, 38, 72, 48, 40… +## $ income <ord> $25000 or more, $20000 - 24999, $25000 or more, $25000 or more… +## $ class <fct> middle class, working class, working class, working class, mid… +## $ finrela <fct> below average, below average, below average, above average, ab… +## $ weight <dbl> 0.8960, 1.0825, 0.5501, 1.0864, 1.0825, 1.0864, 1.0627, 0.4785…

Each row is an individual survey response, containing some basic demographic information on the respondent as well as some additional variables. See ?gss for more information on the variables included and their source. Note that this data (and our examples on it) are for demonstration purposes only, and will not necessarily provide accurate estimates unless weighted properly. For these examples, let’s suppose that this dataset is a representative sample of a population we want to learn about: American adults.

-specify(): Specifying Response (and Explanatory) Variables

+specify(): Specifying Response (and Explanatory) Variables

The specify function can be used to specify which of the variables in the dataset you’re interested in. If you’re only interested in, say, the age of the respondents, you might write:

-
gss %>%
-  specify(response = age)
+
+gss %>%
+  specify(response = age)
## Response: age (numeric)
-## # A tibble: 500 x 1
+## # A tibble: 500 × 1
 ##      age
 ##    <dbl>
 ##  1    36
@@ -177,18 +183,20 @@ 

## 10 33 ## # … with 490 more rows

On the front-end, the output of specify just looks like it selects off the columns in the dataframe that you’ve specified. Checking the class of this object, though:

-
gss %>%
-  specify(response = age) %>%
-  class()
+
+gss %>%
+  specify(response = age) %>%
+  class()
## [1] "infer"      "tbl_df"     "tbl"        "data.frame"

We can see that the infer class has been appended on top of the dataframe classes–this new class stores some extra metadata.

If you’re interested in two variables–age and partyid, for example–you can specify their relationship in one of two (equivalent) ways:

-
# with the named arguments
-gss %>%
-  specify(age ~ partyid)
+
+# as a formula
+gss %>%
+  specify(age ~ partyid)
## Response: age (numeric)
 ## Explanatory: partyid (factor)
-## # A tibble: 500 x 2
+## # A tibble: 500 × 2
 ##      age partyid
 ##    <dbl> <fct>  
 ##  1    36 ind    
@@ -202,12 +210,13 @@ 

## 9 30 rep ## 10 33 dem ## # … with 490 more rows

-
# as a formula
-gss %>%
-  specify(response = age, explanatory = partyid)
+
+# with the named arguments
+gss %>%
+  specify(response = age, explanatory = partyid)
## Response: age (numeric)
 ## Explanatory: partyid (factor)
-## # A tibble: 500 x 2
+## # A tibble: 500 × 2
 ##      age partyid
 ##    <dbl> <fct>  
 ##  1    36 ind    
@@ -222,11 +231,12 @@ 

## 10 33 dem ## # … with 490 more rows

If you’re doing inference on one proportion or a difference in proportions, you will need to use the success argument to specify which level of your response variable is a success. For instance, if you’re interested in the proportion of the population with a college degree, you might use the following code:

-
# specifying for inference on proportions
-gss %>%
-  specify(response = college, success = "degree")
+
+# specifying for inference on proportions
+gss %>%
+  specify(response = college, success = "degree")
## Response: college (factor)
-## # A tibble: 500 x 1
+## # A tibble: 500 × 1
 ##    college  
 ##    <fct>    
 ##  1 degree   
@@ -243,15 +253,16 @@ 

-hypothesize(): Declaring the Null Hypothesis

+hypothesize(): Declaring the Null Hypothesis

The next step in the infer pipeline is often to declare a null hypothesis using hypothesize(). The first step is to supply one of “independence” or “point” to the null argument. If your null hypothesis assumes independence between two variables, then this is all you need to supply to hypothesize():

-
gss %>%
-  specify(college ~ partyid, success = "degree") %>%
-  hypothesize(null = "independence")
+
+gss %>%
+  specify(college ~ partyid, success = "degree") %>%
+  hypothesize(null = "independence")
## Response: college (factor)
 ## Explanatory: partyid (factor)
 ## Null Hypothesis: independence
-## # A tibble: 500 x 2
+## # A tibble: 500 × 2
 ##    college   partyid
 ##    <fct>     <fct>  
 ##  1 degree    ind    
@@ -266,12 +277,13 @@ 

## 10 no degree dem ## # … with 490 more rows

If you’re doing inference on a point estimate, you will also need to provide one of p (the true proportion of successes, between 0 and 1), mu (the true mean), med (the true median), or sigma (the true standard deviation). For instance, if the null hypothesis is that the mean number of hours worked per week in our population is 40, we would write:

-
gss %>%
-  specify(response = hours) %>%
-  hypothesize(null = "point", mu = 40)
+
+gss %>%
+  specify(response = hours) %>%
+  hypothesize(null = "point", mu = 40)
## Response: hours (numeric)
 ## Null Hypothesis: point
-## # A tibble: 500 x 1
+## # A tibble: 500 × 1
 ##    hours
 ##    <dbl>
 ##  1    50
@@ -289,7 +301,7 @@ 

-generate(): Generating the Null Distribution

+generate(): Generating the Null Distribution

Once we’ve asserted our null hypothesis using hypothesize(), we can construct a null distribution based on this hypothesis. We can do this using one of several methods, supplied in the type argument:

  • @@ -302,177 +314,285 @@

    simulate: A value will be sampled from a theoretical distribution with parameters specified in hypothesize() for each replicate. (This option is currently only applicable for testing point estimates.)

Continuing on with our example above, about the average number of hours worked a week, we might write:

-
gss %>%
-  specify(response = hours) %>%
-  hypothesize(null = "point", mu = 40) %>%
-  generate(reps = 1000, type = "bootstrap")
+
+gss %>%
+  specify(response = hours) %>%
+  hypothesize(null = "point", mu = 40) %>%
+  generate(reps = 1000, type = "bootstrap")
## Response: hours (numeric)
 ## Null Hypothesis: point
-## # A tibble: 500,000 x 2
+## # A tibble: 500,000 × 2
 ## # Groups:   replicate [1,000]
 ##    replicate hours
 ##        <int> <dbl>
-##  1         1  48.6
-##  2         1  36.6
-##  3         1  38.6
-##  4         1  73.6
-##  5         1  48.6
-##  6         1  38.6
-##  7         1  61.6
-##  8         1  58.6
-##  9         1  38.6
-## 10         1  23.6
+##  1         1  28.6
+##  2         1  41.6
+##  3         1  54.6
+##  4         1  43.6
+##  5         1  44.6
+##  6         1  43.6
+##  7         1  28.6
+##  8         1  25.6
+##  9         1  78.6
+## 10         1  21.6
 ## # … with 499,990 more rows

In the above example, we take 1000 bootstrap samples to form our null distribution.

To generate a null distribution for the independence of two variables, we could also randomly reshuffle the pairings of explanatory and response variables to break any existing association. For instance, to generate 1000 replicates that can be used to create a null distribution under the assumption that political party affiliation is not affected by age:

-
gss %>%
-  specify(partyid ~ age) %>%
-  hypothesize(null = "independence") %>%
-  generate(reps = 1000, type = "permute")
+
+gss %>%
+  specify(partyid ~ age) %>%
+  hypothesize(null = "independence") %>%
+  generate(reps = 1000, type = "permute")
## Response: partyid (factor)
 ## Explanatory: age (numeric)
 ## Null Hypothesis: independence
-## # A tibble: 500,000 x 3
+## # A tibble: 500,000 × 3
 ## # Groups:   replicate [1,000]
 ##    partyid   age replicate
 ##    <fct>   <dbl>     <int>
-##  1 rep        36         1
-##  2 ind        34         1
+##  1 ind        36         1
+##  2 rep        34         1
 ##  3 dem        24         1
-##  4 ind        42         1
-##  5 rep        31         1
-##  6 rep        32         1
-##  7 dem        48         1
-##  8 ind        36         1
-##  9 rep        30         1
-## 10 dem        33         1
+##  4 rep        42         1
+##  5 dem        31         1
+##  6 dem        32         1
+##  7 ind        48         1
+##  8 rep        36         1
+##  9 ind        30         1
+## 10 ind        33         1
 ## # … with 499,990 more rows

-calculate(): Calculating Summary Statistics

-

Depending on whether you’re carrying out computation-based inference or theory-based inference, you will either supply calculate() with the output of generate() or hypothesize, respectively. The function, for one, takes in a stat argument, which is currently one of “mean”, “median”, “sum”, “sd”, “prop”, “count”, “diff in means”, “diff in medians”, “diff in props”, “Chisq”, “F”, “t”, “z”, “slope”, or “correlation”. For example, continuing our example above to calculate the null distribution of mean hours worked per week:

-
gss %>%
-  specify(response = hours) %>%
-  hypothesize(null = "point", mu = 40) %>%
-  generate(reps = 1000, type = "bootstrap") %>%
-  calculate(stat = "mean")
-
## # A tibble: 1,000 x 2
+calculate(): Calculating Summary Statistics
+

calculate() calculates summary statistics from the output of infer core functions. The function takes in a stat argument, which is currently one of “mean”, “median”, “sum”, “sd”, “prop”, “count”, “diff in means”, “diff in medians”, “diff in props”, “Chisq”, “F”, “t”, “z”, “slope”, or “correlation”. For example, continuing our example above to calculate the null distribution of mean hours worked per week:

+
+gss %>%
+  specify(response = hours) %>%
+  hypothesize(null = "point", mu = 40) %>%
+  generate(reps = 1000, type = "bootstrap") %>%
+  calculate(stat = "mean")
+
## Response: hours (numeric)
+## Null Hypothesis: point
+## # A tibble: 1,000 × 2
 ##    replicate  stat
 ##        <int> <dbl>
-##  1         1  41.1
-##  2         2  38.8
-##  3         3  39.5
-##  4         4  39.5
-##  5         5  39.7
-##  6         6  40.4
-##  7         7  39.6
-##  8         8  41.3
-##  9         9  40.2
-## 10        10  39.8
+##  1         1  39.4
+##  2         2  39.6
+##  3         3  39.2
+##  4         4  40.8
+##  5         5  39.0
+##  6         6  39.4
+##  7         7  39.3
+##  8         8  41.2
+##  9         9  40.5
+## 10        10  40.0
 ## # … with 990 more rows

The output of calculate() here shows us the sample statistic (in this case, the mean) for each of our 1000 replicates. If you’re carrying out inference on differences in means, medians, or proportions, or t and z statistics, you will need to supply an order argument, giving the order in which the explanatory variables should be subtracted. For instance, to find the difference in mean age of those that have a college degree and those that don’t, we might write:

-
gss %>%
-  specify(age ~ college) %>%
-  hypothesize(null = "independence") %>%
-  generate(reps = 1000, type = "permute") %>%
-  calculate("diff in means", order = c("degree", "no degree"))
-
## # A tibble: 1,000 x 2
-##    replicate     stat
-##        <int>    <dbl>
-##  1         1  1.44   
-##  2         2 -0.725  
-##  3         3  0.227  
-##  4         4 -0.00250
-##  5         5 -1.93   
-##  6         6  0.271  
-##  7         7  1.13   
-##  8         8 -1.10   
-##  9         9  1.69   
-## 10        10  1.18   
+
+gss %>%
+  specify(age ~ college) %>%
+  hypothesize(null = "independence") %>%
+  generate(reps = 1000, type = "permute") %>%
+  calculate("diff in means", order = c("degree", "no degree"))
+
## Response: age (numeric)
+## Explanatory: college (factor)
+## Null Hypothesis: independence
+## # A tibble: 1,000 × 2
+##    replicate   stat
+##        <int>  <dbl>
+##  1         1 -0.620
+##  2         2 -0.320
+##  3         3  1.88 
+##  4         4 -1.02 
+##  5         5 -0.531
+##  6         6 -0.972
+##  7         7 -0.478
+##  8         8  1.37 
+##  9         9  1.00 
+## 10        10 -1.55 
 ## # … with 990 more rows

-Other Utilities

-

infer also offers several utilities to extract the meaning out of summary statistics and null distributions—the package provides functions to visualize where a statistic is relative to a distribution (with visualize()), calculate p-values (with get_p_value()), and calculate confidence intervals (with get_confidence_interval()).

+Other Utilities +

infer also offers several utilities to extract the meaning out of summary statistics and distributions—the package provides functions to visualize where a statistic is relative to a distribution (with visualize()), calculate p-values (with get_p_value()), and calculate confidence intervals (with get_confidence_interval()).

To illustrate, we’ll go back to the example of determining whether the mean number of hours worked per week is 40 hours.

-
# find the point estimate
-point_estimate <- gss %>%
-  specify(response = hours) %>%
-  calculate(stat = "mean")
+
+# find the point estimate
+obs_mean <- gss %>%
+  specify(response = hours) %>%
+  calculate(stat = "mean")
 
 # generate a null distribution
-null_dist <- gss %>%
-  specify(response = hours) %>%
-  hypothesize(null = "point", mu = 40) %>%
-  generate(reps = 1000, type = "bootstrap") %>%
-  calculate(stat = "mean")
-

(Notice the warning: Removed 1244 rows containing missing values. This would be worth noting if you were actually carrying out this hypothesis test.)

+null_dist <- gss %>% + specify(response = hours) %>% + hypothesize(null = "point", mu = 40) %>% + generate(reps = 1000, type = "bootstrap") %>% + calculate(stat = "mean")

Our point estimate 41.382 seems pretty close to 40, but a little bit different. We might wonder if this difference is just due to random chance, or if the mean number of hours worked per week in the population really isn’t 40.

We could initially just visualize the null distribution.

-
null_dist %>%
-  visualize()
+
+null_dist %>%
+  visualize()

Where does our sample’s observed statistic lie on this distribution? We can use the obs_stat argument to specify this.

-
null_dist %>%
-  visualize() +
-  shade_p_value(obs_stat = point_estimate, direction = "two-sided")
+
+null_dist %>%
+  visualize() +
+  shade_p_value(obs_stat = obs_mean, direction = "two-sided")

-

Notice that infer has also shaded the regions of the null distribution that are as (or more) extreme than our observed statistic. (Also, note that we now use the + operator to apply the shade_p_value function. This is because visualize outputs a plot object from ggplot2 instead of a data frame, and the + operator is needed to add the p-value layer to the plot object.) The red bar looks like it’s slightly far out on the right tail of the null distribution, so observing a sample mean of 41.382 hours would be somewhat unlikely if the mean was actually 40 hours. How unlikely, though?

-
# get a two-tailed p-value
-p_value <- null_dist %>%
-  get_p_value(obs_stat = point_estimate, direction = "two-sided")
+

Notice that infer has also shaded the regions of the null distribution that are as (or more) extreme than our observed statistic. (Also, note that we now use the + operator to apply the shade_p_value function. This is because visualize outputs a plot object from ggplot2 instead of a data frame, and the + operator is needed to add the p-value layer to the plot object.) The red bar looks like it’s slightly far out on the right tail of the null distribution, so observing a sample mean of 41.382 hours would be somewhat unlikely if the mean was actually 40 hours. How unlikely, though?

+
+# get a two-tailed p-value
+p_value <- null_dist %>%
+  get_p_value(obs_stat = obs_mean, direction = "two-sided")
 
-p_value
-
## # A tibble: 1 x 1
+p_value
+
## # A tibble: 1 × 1
 ##   p_value
 ##     <dbl>
-## 1   0.036
-

It looks like the p-value is 0.036, which is pretty small—if the true mean number of hours worked per week was actually 40, the probability of our sample mean being this far (1.382 hours) from 40 would be 0.036. This may or may not be statistically significantly different, depending on the significance level \(\alpha\) you decided on before you ran this analysis. If you had set \(\alpha = .05\), then this difference would be statistically significant, but if you had set \(\alpha = .01\), then it would not be.

+## 1 0.03
+

It looks like the p-value is 0.03, which is pretty small—if the true mean number of hours worked per week was actually 40, the probability of our sample mean being this far (1.382 hours) from 40 would be 0.03. This may or may not be statistically significantly different, depending on the significance level \(\alpha\) you decided on before you ran this analysis. If you had set \(\alpha = .05\), then this difference would be statistically significant, but if you had set \(\alpha = .01\), then it would not be.

To get a confidence interval around our estimate, we can write:

-
# start with the null distribution
-null_dist %>%
+
+# generate a distribution like the null distribution, 
+# though exclude the null hypothesis from the pipeline
+boot_dist <- gss %>%
+  specify(response = hours) %>%
+  generate(reps = 1000, type = "bootstrap") %>%
+  calculate(stat = "mean")
+
+# start with the bootstrap distribution
+ci <- boot_dist %>%
   # calculate the confidence interval around the point estimate
-  get_confidence_interval(point_estimate = point_estimate,
+  get_confidence_interval(point_estimate = obs_mean,
                           # at the 95% confidence level
-                          level = .95,
+                          level = .95,
                           # using the standard error
-                          type = "se")
-
## # A tibble: 1 x 2
-##   lower upper
-##   <dbl> <dbl>
-## 1  40.1  42.7
-

As you can see, 40 hours per week is not contained in this interval, which aligns with our previous conclusion that this finding is significant at the confidence level \(\alpha = .05\).

+ type = "se") + +ci
+
## # A tibble: 1 × 2
+##   lower_ci upper_ci
+##      <dbl>    <dbl>
+## 1     40.1     42.7
+

As you can see, 40 hours per week is not contained in this interval, which aligns with our previous conclusion that this finding is significant at the confidence level \(\alpha = .05\). To see this interval represented visually, we can use the shade_confidence_interval() utility:

+
+boot_dist %>%
+  visualize() +
+  shade_confidence_interval(endpoints = ci)
+

-Theoretical Methods

-

{infer} also provides functionality to use theoretical methods for "Chisq", "F" and "t" test statistics.

-

Generally, to find a null distribution using theory-based methods, use the same code that you would use to find the null distribution using randomization-based methods, but skip the generate() step. For example, if we wanted to find a null distribution for the relationship between age (age) and party identification (partyid) using randomization, we could write:

-
null_f_distn <- gss %>%
-   specify(age ~ partyid) %>%
-   hypothesize(null = "independence") %>%
-   generate(reps = 1000, type = "permute") %>%
-   calculate(stat = "F")
-

To find the null distribution using theory-based methods, instead, skip the generate() step entirely:

-
null_f_distn_theoretical <- gss %>%
-   specify(age ~ partyid) %>%
-   hypothesize(null = "independence") %>%
-   calculate(stat = "F")
-

We’ll calculate the observed statistic to make use of in the following visualizations—this procedure is the same, regardless of the methods used to find the null distribution.

-
F_hat <- gss %>%
-  specify(age ~ partyid) %>%
-  calculate(stat = "F")
-

Now, instead of just piping the null distribution into visualize(), as we would do if we wanted to visualize the randomization-based null distribution, we also need to provide method = "theoretical" to visualize().

-
visualize(null_f_distn_theoretical, method = "theoretical") +
-  shade_p_value(obs_stat = F_hat, direction = "greater")
-

-

To get a sense of how the theory-based and randomization-based null distributions relate, as well, we can pipe the randomization-based null distribution into visualize() and also specify method = "both"

-
visualize(null_f_distn, method = "both") +
-  shade_p_value(obs_stat = F_hat, direction = "greater")
+Theoretical Methods +

{infer} also provides functionality to use theoretical methods for "Chisq", "F", "t" and "z" distributions.

+

Generally, to find a null distribution using theory-based methods, use the same code that you would use to find the observed statistic elsewhere, replacing calls to calculate() with assume(). For example, to calculate the observed \(t\) statistic (a standardized mean):

+
+# calculate an observed t statistic
+obs_t <- gss %>%
+  specify(response = hours) %>%
+  hypothesize(null = "point", mu = 40) %>%
+  calculate(stat = "t")
+

Then, to define a theoretical \(t\) distribution, we could write:

+
+# switch out calculate with assume to define a distribution
+t_dist <- gss %>%
+  specify(response = hours) %>%
+  assume(distribution = "t")
+

From here, the theoretical distribution interfaces in the same way that simulation-based null distributions do. For example, to interface with p-values:

+
+# visualize the theoretical null distribution
+visualize(t_dist) +
+  shade_p_value(obs_stat = obs_t, direction = "greater")
+

+
+# more exactly, calculate the p-value
+get_p_value(t_dist, obs_t, "greater")
+
## # A tibble: 1 × 1
+##   p_value
+##     <dbl>
+## 1  0.0188
+

Confidence intervals lie on the scale of the data rather than on the standardized scale of the theoretical distribution, so be sure to use the unstandardized observed statistic when working with confidence intervals.

+
+# find the theory-based confidence interval
+theor_ci <- 
+  get_confidence_interval(
+    x = t_dist,
+    level = .95,
+    point_estimate = obs_mean
+  )
+
+theor_ci
+
## # A tibble: 1 × 2
+##   lower_ci upper_ci
+##      <dbl>    <dbl>
+## 1     40.1     42.7
+

When visualized, the \(t\) distribution will be recentered and rescaled to align with the scale of the observed data.

+
+# visualize the theoretical sampling distribution
+visualize(t_dist) +
+  shade_confidence_interval(theor_ci)

-

That’s it! This vignette covers most all of the key functionality of infer. See help(package = "infer") for a full list of functions and vignettes.

+
+
+

+Multiple regression

+

To accommodate randomization-based inference with multiple explanatory variables, the package implements an alternative workflow based on model fitting. Rather than calculate()ing statistics from resampled data, this side of the package allows you to fit() linear models on data resampled according to the null hypothesis, supplying model coefficients for each explanatory variable. For the most part, you can just switch out calculate() for fit() in your calculate()-based workflows.

+

As an example, suppose that we want to fit hours worked per week using the respondent age and college completion status. We could first begin by fitting a linear model to the observed data.

+
+observed_fit <- gss %>%
+  specify(hours ~ age + college) %>%
+  fit()
+

Now, to generate null distributions for each of these terms, we can fit 1000 models to resamples of the gss dataset, where the response hours is permuted in each. Note that this code is the same as the above except for the addition of the hypothesize and generate step.

+
+null_fits <- gss %>%
+  specify(hours ~ age + college) %>%
+  hypothesize(null = "independence") %>%
+  generate(reps = 1000, type = "permute") %>%
+  fit()
+
+null_fits
+
## # A tibble: 3,000 × 3
+## # Groups:   replicate [1,000]
+##    replicate term          estimate
+##        <int> <chr>            <dbl>
+##  1         1 intercept     41.1    
+##  2         1 age            0.00693
+##  3         1 collegedegree  0.0422 
+##  4         2 intercept     44.2    
+##  5         2 age           -0.0773 
+##  6         2 collegedegree  0.924  
+##  7         3 intercept     42.8    
+##  8         3 age           -0.0370 
+##  9         3 collegedegree  0.260  
+## 10         4 intercept     44.2    
+## # … with 2,990 more rows
+

To permute variables other than the response variable, the variables argument to generate() allows you to choose columns from the data to permute. Note that any derived effects that depend on these columns (e.g., interaction effects) will also be affected.

+

Beyond this point, observed fits and distributions from null fits interface exactly like analogous outputs from calculate(). For instance, we can use the following code to calculate a 95% confidence interval from these objects.

+
+get_confidence_interval(
+  null_fits, 
+  point_estimate = observed_fit, 
+  level = .95
+)
+
## # A tibble: 3 × 3
+##   term          lower_ci upper_ci
+##   <chr>            <dbl>    <dbl>
+## 1 age            -0.0950   0.0939
+## 2 collegedegree  -2.62     2.62  
+## 3 intercept      37.5     45.2
+

Or, we can shade p-values for each of these observed regression coefficients from the observed data.

+
+visualize(null_fits) + 
+  shade_p_value(observed_fit, direction = "both")
+

+
+
+

+Conclusion

+

That’s it! This vignette covers most all of the key functionality of infer. See help(package = "infer") for a full list of functions and vignettes.

@@ -484,32 +604,24 @@

-
- - + + + + + diff --git a/docs/articles/infer_files/figure-html/unnamed-chunk-10-1.png b/docs/articles/infer_files/figure-html/unnamed-chunk-10-1.png new file mode 100644 index 00000000..79bd4b74 Binary files /dev/null and b/docs/articles/infer_files/figure-html/unnamed-chunk-10-1.png differ diff --git a/docs/articles/infer_files/figure-html/unnamed-chunk-4-1.png b/docs/articles/infer_files/figure-html/unnamed-chunk-4-1.png new file mode 100644 index 00000000..9a05e6ea Binary files /dev/null and b/docs/articles/infer_files/figure-html/unnamed-chunk-4-1.png differ diff --git a/docs/articles/infer_files/figure-html/unnamed-chunk-6-1.png b/docs/articles/infer_files/figure-html/unnamed-chunk-6-1.png index 0b125c0e..98c6c7cd 100644 Binary files a/docs/articles/infer_files/figure-html/unnamed-chunk-6-1.png and b/docs/articles/infer_files/figure-html/unnamed-chunk-6-1.png differ diff --git a/docs/articles/infer_files/figure-html/visualize-1.png b/docs/articles/infer_files/figure-html/visualize-1.png index 7d30e8f6..bcb7a261 100644 Binary files a/docs/articles/infer_files/figure-html/visualize-1.png and b/docs/articles/infer_files/figure-html/visualize-1.png differ diff --git a/docs/articles/infer_files/figure-html/visualize-ci-1.png b/docs/articles/infer_files/figure-html/visualize-ci-1.png new file mode 100644 index 00000000..22375c58 Binary files /dev/null and b/docs/articles/infer_files/figure-html/visualize-ci-1.png differ diff --git a/docs/articles/infer_files/figure-html/visualize2-1.png b/docs/articles/infer_files/figure-html/visualize2-1.png index be5c95f5..bb70bfc8 100644 Binary files a/docs/articles/infer_files/figure-html/visualize2-1.png and b/docs/articles/infer_files/figure-html/visualize2-1.png differ diff --git a/docs/articles/infer_files/header-attrs-2.8/header-attrs.js b/docs/articles/infer_files/header-attrs-2.8/header-attrs.js new file mode 100644 index 00000000..dd57d92e --- /dev/null +++ b/docs/articles/infer_files/header-attrs-2.8/header-attrs.js @@ -0,0 +1,12 @@ +// Pandoc 2.9 adds attributes on both header and div. We remove the former (to +// be compatible with the behavior of Pandoc < 2.8). +document.addEventListener('DOMContentLoaded', function(e) { + var hs = document.querySelectorAll("div.section[class*='level'] > :first-child"); + var i, h, a; + for (i = 0; i < hs.length; i++) { + h = hs[i]; + if (!/^h[1-6]$/i.test(h.tagName)) continue; // it should be a header h1-h6 + a = h.attributes; + while (a.length > 0) h.removeAttribute(a[0].name); + } +}); diff --git a/docs/articles/observed_stat_examples.html b/docs/articles/observed_stat_examples.html index 7fc7abf1..6f613e3a 100644 --- a/docs/articles/observed_stat_examples.html +++ b/docs/articles/observed_stat_examples.html @@ -12,21 +12,30 @@ - + - - - - + + + + + + + +
@@ -105,13 +109,13 @@ -
+
@@ -120,856 +124,1269 @@

Full infer Pipeline Examples

-Introduction

+Introduction

This vignette is intended to provide a set of examples that nearly exhaustively demonstrate the functionalities provided by infer. Commentary on these examples is limited—for more discussion of the intuition behind the package, see the “Getting to Know infer” vignette, accessible by calling vignette("infer").

Throughout this vignette, we’ll make use of the gss dataset supplied by infer, which contains a sample of data from the General Social Survey. See ?gss for more information on the variables included and their source. Note that this data (and our examples on it) are for demonstration purposes only, and will not necessarily provide accurate estimates unless weighted properly. For these examples, let’s suppose that this dataset is a representative sample of a population we want to learn about: American adults. The data looks like this:

-
# load in the dataset
-data(gss)
+
+# load in the dataset
+data(gss)
 
 # take a look at its structure
-dplyr::glimpse(gss)
+dplyr::glimpse(gss)
## Rows: 500
 ## Columns: 11
-## $ year    <dbl> 2014, 1994, 1998, 1996, 1994, 1996, 1990, 2016, 2000, 1998, 2…
-## $ age     <dbl> 36, 34, 24, 42, 31, 32, 48, 36, 30, 33, 21, 30, 38, 49, 25, 5…
-## $ sex     <fct> male, female, male, male, male, female, female, female, femal…
-## $ college <fct> degree, no degree, degree, no degree, degree, no degree, no d…
-## $ partyid <fct> ind, rep, ind, ind, rep, rep, dem, ind, rep, dem, dem, ind, d…
-## $ hompop  <dbl> 3, 4, 1, 4, 2, 4, 2, 1, 5, 2, 4, 3, 4, 4, 2, 2, 3, 2, 1, 2, 5…
-## $ hours   <dbl> 50, 31, 40, 40, 40, 53, 32, 20, 40, 40, 23, 52, 38, 72, 48, 4…
-## $ income  <ord> $25000 or more, $20000 - 24999, $25000 or more, $25000 or mor…
-## $ class   <fct> middle class, working class, working class, working class, mi…
-## $ finrela <fct> below average, below average, below average, above average, a…
-## $ weight  <dbl> 0.8960, 1.0825, 0.5501, 1.0864, 1.0825, 1.0864, 1.0627, 0.478…
+## $ year <dbl> 2014, 1994, 1998, 1996, 1994, 1996, 1990, 2016, 2000, 1998, 20… +## $ age <dbl> 36, 34, 24, 42, 31, 32, 48, 36, 30, 33, 21, 30, 38, 49, 25, 56… +## $ sex <fct> male, female, male, male, male, female, female, female, female… +## $ college <fct> degree, no degree, degree, no degree, degree, no degree, no de… +## $ partyid <fct> ind, rep, ind, ind, rep, rep, dem, ind, rep, dem, dem, ind, de… +## $ hompop <dbl> 3, 4, 1, 4, 2, 4, 2, 1, 5, 2, 4, 3, 4, 4, 2, 2, 3, 2, 1, 2, 5,… +## $ hours <dbl> 50, 31, 40, 40, 40, 53, 32, 20, 40, 40, 23, 52, 38, 72, 48, 40… +## $ income <ord> $25000 or more, $20000 - 24999, $25000 or more, $25000 or more… +## $ class <fct> middle class, working class, working class, working class, mid… +## $ finrela <fct> below average, below average, below average, above average, ab… +## $ weight <dbl> 0.8960, 1.0825, 0.5501, 1.0864, 1.0825, 1.0864, 1.0627, 0.4785…

-Hypothesis tests

+Hypothesis tests

-One numerical variable (mean)

+One numerical variable (mean)

Calculating the observed statistic,

-
x_bar <- gss %>%
-  specify(response = hours) %>%
-  calculate(stat = "mean")
+
+x_bar <- gss %>%
+  specify(response = hours) %>%
+  calculate(stat = "mean")
+

Alternatively, using the observe() wrapper to calculate the observed statistic,

+
+x_bar <- gss %>%
+  observe(response = hours, stat = "mean")

Then, generating the null distribution,

-
null_distn <- gss %>%
-  specify(response = hours) %>%
-  hypothesize(null = "point", mu = 40) %>%
-  generate(reps = 1000) %>%
-  calculate(stat = "mean")
+
+null_dist <- gss %>%
+  specify(response = hours) %>%
+  hypothesize(null = "point", mu = 40) %>%
+  generate(reps = 1000) %>%
+  calculate(stat = "mean")

Visualizing the observed statistic alongside the null distribution,

-
visualize(null_distn) +
-  shade_p_value(obs_stat = x_bar, direction = "two-sided")
-

+
+visualize(null_dist) +
+  shade_p_value(obs_stat = x_bar, direction = "two-sided")
+

Calculating the p-value from the null distribution and observed statistic,

-
null_distn %>%
-  get_p_value(obs_stat = x_bar, direction = "two-sided")
-
## # A tibble: 1 x 1
+
+null_dist %>%
+  get_p_value(obs_stat = x_bar, direction = "two-sided")
+
## # A tibble: 1 × 1
 ##   p_value
 ##     <dbl>
-## 1   0.032
+## 1 0.026

-One numerical variable (standardized mean \(t\))

+One numerical variable (standardized mean \(t\))

Calculating the observed statistic,

-
t_bar <- gss %>%
-  specify(response = hours) %>%
-  hypothesize(null = "point", mu = 40) %>%
-  calculate(stat = "t")
-

Alternatively, using the wrapper to calculate the test statistic,

-
t_bar <- gss %>%
-  t_stat(response = hours, mu = 40)
+
+t_bar <- gss %>%
+  specify(response = hours) %>%
+  hypothesize(null = "point", mu = 40) %>%
+  calculate(stat = "t")
+

Alternatively, using the observe() wrapper to calculate the observed statistic,

+
+t_bar <- gss %>%
+  observe(response = hours, null = "point", mu = 40, stat = "t")

Then, generating the null distribution,

-
null_distn <- gss %>%
-  specify(response = hours) %>%
-  hypothesize(null = "point", mu = 40) %>%
-  generate(reps = 1000) %>%
-  calculate(stat = "t")
-

Alternatively, finding the null distribution using theoretical methods by skipping the generate() step,

-
null_distn_theoretical <- gss %>%
-  specify(response = hours) %>%
-  hypothesize(null = "point", mu = 40) %>%
-  calculate(stat = "t")
+
+null_dist <- gss %>%
+  specify(response = hours) %>%
+  hypothesize(null = "point", mu = 40) %>%
+  generate(reps = 1000) %>%
+  calculate(stat = "t")
+

Alternatively, finding the null distribution using theoretical methods using the assume() verb,

+
+null_dist_theory <- gss %>%
+  specify(response = hours)  %>%
+  assume("t")

Visualizing the observed statistic alongside the null distribution,

-
visualize(null_distn) +
-  shade_p_value(obs_stat = t_bar, direction = "two-sided")
-

-

Alternatively, visualizing the observed statistic using the theory-based null distribution,

-
visualize(null_distn_theoretical, method = "theoretical") +
-  shade_p_value(obs_stat = t_bar, direction = "two-sided")
+
+visualize(null_dist) +
+  shade_p_value(obs_stat = t_bar, direction = "two-sided")

-

Alternatively, visualizing the observed statistic using both of the null distributions,

-
visualize(null_distn, method = "both") +
-  shade_p_value(obs_stat = t_bar, direction = "two-sided")
+

Alternatively, visualizing the observed statistic using the theory-based null distribution,

+
+visualize(null_dist_theory) +
+  shade_p_value(obs_stat = t_bar, direction = "two-sided")

+

Alternatively, visualizing the observed statistic using both of the null distributions,

+
+visualize(null_dist, method = "both") +
+  shade_p_value(obs_stat = t_bar, direction = "two-sided")
+

Note that the above code makes use of the randomization-based null distribution.

Calculating the p-value from the null distribution and observed statistic,

-
null_distn %>%
-  get_p_value(obs_stat = t_bar, direction = "two-sided")
-
## # A tibble: 1 x 1
+
+null_dist %>%
+  get_p_value(obs_stat = t_bar, direction = "two-sided")
+
## # A tibble: 1 × 1
 ##   p_value
 ##     <dbl>
 ## 1   0.044

Alternatively, using the t_test wrapper:

-
gss %>%
-  t_test(response = hours, mu = 40)
-
## # A tibble: 1 x 6
-##   statistic  t_df p_value alternative lower_ci upper_ci
-##       <dbl> <dbl>   <dbl> <chr>          <dbl>    <dbl>
-## 1      2.09   499  0.0376 two.sided       40.1     42.7
+
+gss %>%
+  t_test(response = hours, mu = 40)
+
## # A tibble: 1 × 7
+##   statistic  t_df p_value alternative estimate lower_ci upper_ci
+##       <dbl> <dbl>   <dbl> <chr>          <dbl>    <dbl>    <dbl>
+## 1      2.09   499  0.0376 two.sided       41.4     40.1     42.7
+

infer does not support testing on one numerical variable via the z distribution.

-One numerical variable (median)

+One numerical variable (median)

Calculating the observed statistic,

-
x_tilde <- gss %>%
-  specify(response = age) %>%
-  calculate(stat = "median")
+
+x_tilde <- gss %>%
+  specify(response = age) %>%
+  calculate(stat = "median")
+

Alternatively, using the observe() wrapper to calculate the observed statistic,

+
+x_tilde <- gss %>%
+  observe(response = age, stat = "median")

Then, generating the null distribution,

-
null_distn <- gss %>%
-  specify(response = age) %>%
-  hypothesize(null = "point", med = 40) %>%
-  generate(reps = 1000) %>%
-  calculate(stat = "median")
+
+null_dist <- gss %>%
+  specify(response = age) %>%
+  hypothesize(null = "point", med = 40) %>% 
+  generate(reps = 1000) %>% 
+  calculate(stat = "median")

Visualizing the observed statistic alongside the null distribution,

-
visualize(null_distn) +
-  shade_p_value(obs_stat = x_tilde, direction = "two-sided")
-

+
+visualize(null_dist) +
+  shade_p_value(obs_stat = x_tilde, direction = "two-sided")
+

Calculating the p-value from the null distribution and observed statistic,

-
null_distn %>%
-  get_p_value(obs_stat = x_tilde, direction = "two-sided")
-
## # A tibble: 1 x 1
+
+null_dist %>%
+  get_p_value(obs_stat = x_tilde, direction = "two-sided")
+
## # A tibble: 1 × 1
 ##   p_value
 ##     <dbl>
 ## 1   0.008

-One categorical (one proportion)

+One categorical (one proportion)

Calculating the observed statistic,

-
p_hat <- gss %>%
-  specify(response = sex, success = "female") %>%
-  calculate(stat = "prop")
+
+p_hat <- gss %>%
+  specify(response = sex, success = "female") %>%
+  calculate(stat = "prop")
+

Alternatively, using the observe() wrapper to calculate the observed statistic,

+
+p_hat <- gss %>%
+  observe(response = sex, success = "female", stat = "prop")

Then, generating the null distribution,

-
null_distn <- gss %>%
-  specify(response = sex, success = "female") %>%
-  hypothesize(null = "point", p = .5) %>%
-  generate(reps = 1000) %>%
-  calculate(stat = "prop")
+
+null_dist <- gss %>%
+  specify(response = sex, success = "female") %>%
+  hypothesize(null = "point", p = .5) %>%
+  generate(reps = 1000) %>%
+  calculate(stat = "prop")

Visualizing the observed statistic alongside the null distribution,

-
visualize(null_distn) +
-  shade_p_value(obs_stat = p_hat, direction = "two-sided")
-

+
+visualize(null_dist) +
+  shade_p_value(obs_stat = p_hat, direction = "two-sided")
+

Calculating the p-value from the null distribution and observed statistic,

-
null_distn %>%
-  get_p_value(obs_stat = p_hat, direction = "two-sided")
-
## # A tibble: 1 x 1
+
+null_dist %>%
+  get_p_value(obs_stat = p_hat, direction = "two-sided")
+
## # A tibble: 1 × 1
 ##   p_value
 ##     <dbl>
-## 1   0.274
+## 1 0.248

Note that logical variables will be coerced to factors:

-
null_distn <- gss %>%
-  dplyr::mutate(is_female = (sex == "female")) %>%
-  specify(response = is_female, success = "TRUE") %>%
-  hypothesize(null = "point", p = .5) %>%
-  generate(reps = 1000) %>%
-  calculate(stat = "prop")
+
+null_dist <- gss %>%
+  dplyr::mutate(is_female = (sex == "female")) %>%
+  specify(response = is_female, success = "TRUE") %>%
+  hypothesize(null = "point", p = .5) %>%
+  generate(reps = 1000) %>%
+  calculate(stat = "prop")

-One categorical variable (standardized proportion \(z\))

-

While the standardized proportion \(z\) statistic has not yet been implemented in the randomization-based framework, the package supplies a wrapper around prop.test to allow for tests of a single proportion on tidy data.

-
prop_test(gss,
-          college ~ NULL,
-          p = .2)
-
## # A tibble: 1 x 4
+One categorical variable (standardized proportion \(z\))
+

Calculating the observed statistic,

+
+p_hat <- gss %>%
+  specify(response = sex, success = "female") %>%
+  hypothesize(null = "point", p = .5) %>%
+  calculate(stat = "z")
+

Alternatively, using the observe() wrapper to calculate the observed statistic,

+
+p_hat <- gss %>%
+  observe(response = sex, success = "female", null = "point", p = .5, stat = "z")
+

Then, generating the null distribution,

+
+null_dist <- gss %>%
+  specify(response = sex, success = "female") %>%
+  hypothesize(null = "point", p = .5) %>%
+  generate(reps = 1000, type = "draw") %>%
+  calculate(stat = "z")
+

Visualizing the observed statistic alongside the null distribution,

+
+visualize(null_dist) +
+  shade_p_value(obs_stat = p_hat, direction = "two-sided")
+

+

Calculating the p-value from the null distribution and observed statistic,

+
+null_dist %>%
+  get_p_value(obs_stat = p_hat, direction = "two-sided")
+
## # A tibble: 1 × 1
+##   p_value
+##     <dbl>
+## 1   0.244
+

The package also supplies a wrapper around prop.test for tests of a single proportion on tidy data.

+
+prop_test(gss,
+          college ~ NULL,
+          p = .2)
+
## # A tibble: 1 × 4
 ##   statistic chisq_df   p_value alternative
 ##       <dbl>    <int>     <dbl> <chr>      
 ## 1      636.        1 2.98e-140 two.sided
+

infer does not support testing two means via the z distribution.

-Two categorical (2 level) variables

+Two categorical (2 level) variables

The infer package provides several statistics to work with data of this type. One of them is the statistic for difference in proportions.

Calculating the observed statistic,

-
d_hat <- gss %>%
-  specify(college ~ sex, success = "no degree") %>%
-  calculate(stat = "diff in props", order = c("female", "male"))
+
+d_hat <- gss %>% 
+  specify(college ~ sex, success = "no degree") %>%
+  calculate(stat = "diff in props", order = c("female", "male"))
+

Alternatively, using the observe() wrapper to calculate the observed statistic,

+
+d_hat <- gss %>% 
+  observe(college ~ sex, success = "no degree", 
+          stat = "diff in props", order = c("female", "male"))

Then, generating the null distribution,

-
null_distn <- gss %>%
-  specify(college ~ sex, success = "no degree") %>%
-  hypothesize(null = "independence") %>%
-  generate(reps = 1000) %>%
-  calculate(stat = "diff in props", order = c("female", "male"))
+
+null_dist <- gss %>%
+  specify(college ~ sex, success = "no degree") %>%
+  hypothesize(null = "independence") %>% 
+  generate(reps = 1000) %>% 
+  calculate(stat = "diff in props", order = c("female", "male"))

Visualizing the observed statistic alongside the null distribution,

-
visualize(null_distn) +
-  shade_p_value(obs_stat = d_hat, direction = "two-sided")
-

+
+visualize(null_dist) +
+  shade_p_value(obs_stat = d_hat, direction = "two-sided")
+

Calculating the p-value from the null distribution and observed statistic,

-
null_distn %>%
-  get_p_value(obs_stat = d_hat, direction = "two-sided")
-
## # A tibble: 1 x 1
+
+null_dist %>%
+  get_p_value(obs_stat = d_hat, direction = "two-sided")
+
## # A tibble: 1 × 1
 ##   p_value
 ##     <dbl>
 ## 1       1

infer also provides functionality to calculate ratios of proportions. The workflow looks similar to that for diff in props.

Calculating the observed statistic,

-
r_hat <- gss %>%
-  specify(college ~ sex, success = "no degree") %>%
-  calculate(stat = "ratio of props", order = c("female", "male"))
+
+r_hat <- gss %>% 
+  specify(college ~ sex, success = "no degree") %>%
+  calculate(stat = "ratio of props", order = c("female", "male"))
+

Alternatively, using the observe() wrapper to calculate the observed statistic,

+
+r_hat <- gss %>% 
+  observe(college ~ sex, success = "no degree",
+          stat = "ratio of props", order = c("female", "male"))

Then, generating the null distribution,

-
null_distn <- gss %>%
-  specify(college ~ sex, success = "no degree") %>%
-  hypothesize(null = "independence") %>%
-  generate(reps = 1000) %>%
-  calculate(stat = "ratio of props", order = c("female", "male"))
+
+null_dist <- gss %>%
+  specify(college ~ sex, success = "no degree") %>%
+  hypothesize(null = "independence") %>% 
+  generate(reps = 1000) %>% 
+  calculate(stat = "ratio of props", order = c("female", "male"))

Visualizing the observed statistic alongside the null distribution,

-
visualize(null_distn) +
-  shade_p_value(obs_stat = r_hat, direction = "two-sided")
-

+
+visualize(null_dist) +
+  shade_p_value(obs_stat = r_hat, direction = "two-sided")
+

Calculating the p-value from the null distribution and observed statistic,

-
null_distn %>%
-  get_p_value(obs_stat = r_hat, direction = "two-sided")
-
## # A tibble: 1 x 1
+
+null_dist %>%
+  get_p_value(obs_stat = r_hat, direction = "two-sided")
+
## # A tibble: 1 × 1
 ##   p_value
 ##     <dbl>
-## 1   0.964
+## 1 1

In addition, the package provides functionality to calculate odds ratios. The workflow also looks similar to that for diff in props.

Calculating the observed statistic,

-
or_hat <- gss %>%
-  specify(college ~ sex, success = "no degree") %>%
-  calculate(stat = "odds ratio", order = c("female", "male"))
+
+or_hat <- gss %>% 
+  specify(college ~ sex, success = "no degree") %>%
+  calculate(stat = "odds ratio", order = c("female", "male"))

Then, generating the null distribution,

-
null_distn <- gss %>%
-  specify(college ~ sex, success = "no degree") %>%
-  hypothesize(null = "independence") %>%
-  generate(reps = 1000) %>%
-  calculate(stat = "odds ratio", order = c("female", "male"))
+
+null_dist <- gss %>%
+  specify(college ~ sex, success = "no degree") %>%
+  hypothesize(null = "independence") %>% 
+  generate(reps = 1000) %>% 
+  calculate(stat = "odds ratio", order = c("female", "male"))

Visualizing the observed statistic alongside the null distribution,

-
visualize(null_distn) +
-  shade_p_value(obs_stat = or_hat, direction = "two-sided")
-

+
+visualize(null_dist) +
+  shade_p_value(obs_stat = or_hat, direction = "two-sided")
+

Calculating the p-value from the null distribution and observed statistic,

-
null_distn %>%
-  get_p_value(obs_stat = or_hat, direction = "two-sided")
-
## # A tibble: 1 x 1
+
+null_dist %>%
+  get_p_value(obs_stat = or_hat, direction = "two-sided")
+
## # A tibble: 1 × 1
 ##   p_value
 ##     <dbl>
 ## 1       1

-Two categorical (2 level) variables (z)

+Two categorical (2 level) variables (z)

Finding the standardized observed statistic,

-
z_hat <- gss %>%
-  specify(college ~ sex, success = "no degree") %>%
-  calculate(stat = "z", order = c("female", "male"))
+
+z_hat <- gss %>% 
+  specify(college ~ sex, success = "no degree") %>%
+  hypothesize(null = "independence") %>%
+  calculate(stat = "z", order = c("female", "male"))
+

Alternatively, using the observe() wrapper to calculate the observed statistic,

+
+z_hat <- gss %>% 
+  observe(college ~ sex, success = "no degree",
+          stat = "z", order = c("female", "male"))

Then, generating the null distribution,

-
null_distn <- gss %>%
-  specify(college ~ sex, success = "no degree") %>%
-  hypothesize(null = "independence") %>%
-  generate(reps = 1000) %>%
-  calculate(stat = "z", order = c("female", "male"))
-

Alternatively, finding the null distribution using theoretical methods by skipping the generate() step,

-
null_distn_theoretical <- gss %>%
-  specify(college ~ sex, success = "no degree") %>%
-  hypothesize(null = "independence") %>%
-  calculate(stat = "z", order = c("female", "male"))
+
+null_dist <- gss %>%
+  specify(college ~ sex, success = "no degree") %>%
+  hypothesize(null = "independence") %>% 
+  generate(reps = 1000) %>% 
+  calculate(stat = "z", order = c("female", "male"))
+

Alternatively, finding the null distribution using theoretical methods using the assume() verb,

+
+null_dist_theory <- gss %>%
+  specify(college ~ sex, success = "no degree") %>%
+  assume("z")

Visualizing the observed statistic alongside the null distribution,

-
visualize(null_distn) +
-  shade_p_value(obs_stat = z_hat, direction = "two-sided")
-

+
+visualize(null_dist) +
+  shade_p_value(obs_stat = z_hat, direction = "two-sided")
+

Alternatively, visualizing the observed statistic using the theory-based null distribution,

-
visualize(null_distn_theoretical, method = "theoretical") +
-  shade_p_value(obs_stat = z_hat, direction = "two-sided")
-

+
+visualize(null_dist_theory) +
+  shade_p_value(obs_stat = z_hat, direction = "two-sided")
+

Alternatively, visualizing the observed statistic using both of the null distributions,

-
visualize(null_distn, method = "both") +
-  shade_p_value(obs_stat = z_hat, direction = "two-sided")
-

+
+visualize(null_dist, method = "both") +
+  shade_p_value(obs_stat = z_hat, direction = "two-sided")
+

Note that the above code makes use of the randomization-based null distribution.

Calculating the p-value from the null distribution and observed statistic,

-
null_distn %>%
-  get_p_value(obs_stat = z_hat, direction = "two-sided")
-
## # A tibble: 1 x 1
+
+null_dist %>%
+  get_p_value(obs_stat = z_hat, direction = "two-sided")
+
## # A tibble: 1 × 1
 ##   p_value
 ##     <dbl>
 ## 1       1

Note the similarities in this plot and the previous one.

The package also supplies a wrapper around prop.test to allow for tests of equality of proportions on tidy data.

-
prop_test(gss,
-          college ~ sex,
-          order = c("female", "male"))
-
## # A tibble: 1 x 6
+
+prop_test(gss, 
+          college ~ sex,  
+          order = c("female", "male"))
+
## # A tibble: 1 × 6
 ##   statistic chisq_df p_value alternative lower_ci upper_ci
 ##       <dbl>    <dbl>   <dbl> <chr>          <dbl>    <dbl>
 ## 1 0.0000204        1   0.996 two.sided     -0.101   0.0917

-One categorical (>2 level) - GoF

+One categorical (>2 level) - GoF

Calculating the observed statistic,

Note the need to add in the hypothesized values here to compute the observed statistic.

-
Chisq_hat <- gss %>%
-  specify(response = finrela) %>%
-  hypothesize(null = "point",
-              p = c("far below average" = 1/6,
-                    "below average" = 1/6,
-                    "average" = 1/6,
-                    "above average" = 1/6,
-                    "far above average" = 1/6,
-                    "DK" = 1/6)) %>%
-  calculate(stat = "Chisq")
-

Alternatively, using the chisq_stat wrapper to calculate the test statistic,

-
Chisq_hat <- gss %>%
-  chisq_stat(response = finrela,
-             p = c("far below average" = 1/6,
-                   "below average" = 1/6,
-                   "average" = 1/6,
-                   "above average" = 1/6,
-                   "far above average" = 1/6,
-                   "DK" = 1/6))
+
+Chisq_hat <- gss %>%
+  specify(response = finrela) %>%
+  hypothesize(null = "point",
+              p = c("far below average" = 1/6,
+                    "below average" = 1/6,
+                    "average" = 1/6,
+                    "above average" = 1/6,
+                    "far above average" = 1/6,
+                    "DK" = 1/6)) %>%
+  calculate(stat = "Chisq")
+

Alternatively, using the observe() wrapper to calculate the observed statistic,

+
+Chisq_hat <- gss %>%
+  observe(response = finrela,
+          null = "point",
+          p = c("far below average" = 1/6,
+                "below average" = 1/6,
+                "average" = 1/6,
+                "above average" = 1/6,
+                "far above average" = 1/6,
+                "DK" = 1/6),
+          stat = "Chisq")

Then, generating the null distribution,

-
null_distn <- gss %>%
-  specify(response = finrela) %>%
-  hypothesize(null = "point",
-              p = c("far below average" = 1/6,
-                    "below average" = 1/6,
-                    "average" = 1/6,
-                    "above average" = 1/6,
-                    "far above average" = 1/6,
-                    "DK" = 1/6)) %>%
-  generate(reps = 1000, type = "simulate") %>%
-  calculate(stat = "Chisq")
-

Alternatively, finding the null distribution using theoretical methods by skipping the generate() step,

-
null_distn_theoretical <- gss %>%
-  specify(response = finrela) %>%
-  hypothesize(null = "point",
-              p = c("far below average" = 1/6,
-                    "below average" = 1/6,
-                    "average" = 1/6,
-                    "above average" = 1/6,
-                    "far above average" = 1/6,
-                    "DK" = 1/6)) %>%
-  calculate(stat = "Chisq")
+
+null_dist <- gss %>%
+  specify(response = finrela) %>%
+  hypothesize(null = "point",
+              p = c("far below average" = 1/6,
+                    "below average" = 1/6,
+                    "average" = 1/6,
+                    "above average" = 1/6,
+                    "far above average" = 1/6,
+                    "DK" = 1/6)) %>%
+  generate(reps = 1000, type = "draw") %>%
+  calculate(stat = "Chisq")
+

Alternatively, finding the null distribution using theoretical methods using the assume() verb,

+
+null_dist_theory <- gss %>%
+  specify(response = finrela) %>%
+  assume("Chisq")

Visualizing the observed statistic alongside the null distribution,

-
visualize(null_distn) +
-  shade_p_value(obs_stat = Chisq_hat, direction = "greater")
-

+
+visualize(null_dist) +
+  shade_p_value(obs_stat = Chisq_hat, direction = "greater")
+

Alternatively, visualizing the observed statistic using the theory-based null distribution,

-
visualize(null_distn_theoretical, method = "theoretical") +
-  shade_p_value(obs_stat = Chisq_hat, direction = "greater")
-

+
+visualize(null_dist_theory) +
+  shade_p_value(obs_stat = Chisq_hat, direction = "greater")
+

Alternatively, visualizing the observed statistic using both of the null distributions,

-
visualize(null_distn_theoretical, method = "both") +
-  shade_p_value(obs_stat = Chisq_hat, direction = "greater")
-

+
+visualize(null_dist_theory, method = "both") +
+  shade_p_value(obs_stat = Chisq_hat, direction = "greater")
+

Note that the above code makes use of the randomization-based null distribution.

Calculating the p-value from the null distribution and observed statistic,

-
null_distn %>%
-  get_p_value(obs_stat = Chisq_hat, direction = "greater")
-
## # A tibble: 1 x 1
+
+null_dist %>%
+  get_p_value(obs_stat = Chisq_hat, direction = "greater")
+
## # A tibble: 1 × 1
 ##   p_value
 ##     <dbl>
 ## 1       0

Alternatively, using the chisq_test wrapper:

-
chisq_test(gss,
-           response = finrela,
-           p = c("far below average" = 1/6,
-                 "below average" = 1/6,
-                 "average" = 1/6,
-                 "above average" = 1/6,
-                 "far above average" = 1/6,
-                 "DK" = 1/6))
-
## # A tibble: 1 x 3
+
+chisq_test(gss, 
+           response = finrela,
+           p = c("far below average" = 1/6,
+                 "below average" = 1/6,
+                 "average" = 1/6,
+                 "above average" = 1/6,
+                 "far above average" = 1/6,
+                 "DK" = 1/6))
+
## # A tibble: 1 × 3
 ##   statistic chisq_df   p_value
 ##       <dbl>    <dbl>     <dbl>
 ## 1      488.        5 3.13e-103

-Two categorical (>2 level): Chi-squared test of independence

+Two categorical (>2 level): Chi-squared test of independence

Calculating the observed statistic,

-
Chisq_hat <- gss %>%
-  specify(formula = finrela ~ sex) %>%
-  calculate(stat = "Chisq")
-

Alternatively, using the wrapper to calculate the test statistic,

-
Chisq_hat <- gss %>%
-  chisq_stat(formula = finrela ~ sex)
+
+Chisq_hat <- gss %>%
+  specify(formula = finrela ~ sex) %>% 
+  hypothesize(null = "independence") %>%
+  calculate(stat = "Chisq")
+

Alternatively, using the observe() wrapper to calculate the observed statistic,

+
+Chisq_hat <- gss %>%
+  observe(formula = finrela ~ sex, stat = "Chisq")

Then, generating the null distribution,

-
null_distn <- gss %>%
-  specify(finrela ~ sex) %>%
-  hypothesize(null = "independence") %>%
-  generate(reps = 1000, type = "permute") %>%
-  calculate(stat = "Chisq")
-

Alternatively, finding the null distribution using theoretical methods by skipping the generate() step,

-
null_distn_theoretical <- gss %>%
-  specify(finrela ~ sex) %>%
-  hypothesize(null = "independence") %>%
-  calculate(stat = "Chisq")
+
+null_dist <- gss %>%
+  specify(finrela ~ sex) %>%
+  hypothesize(null = "independence") %>% 
+  generate(reps = 1000, type = "permute") %>% 
+  calculate(stat = "Chisq")
+

Alternatively, finding the null distribution using theoretical methods using the assume() verb,

+
+null_dist_theory <- gss %>%
+  specify(finrela ~ sex) %>%
+  assume(distribution = "Chisq")

Visualizing the observed statistic alongside the null distribution,

-
visualize(null_distn) +
-  shade_p_value(obs_stat = Chisq_hat, direction = "greater")
-

+
+visualize(null_dist) +
+  shade_p_value(obs_stat = Chisq_hat, direction = "greater")
+

Alternatively, visualizing the observed statistic using the theory-based null distribution,

-
visualize(null_distn_theoretical, method = "theoretical") +
-  shade_p_value(obs_stat = Chisq_hat, direction = "greater")
-

+
+visualize(null_dist_theory) +
+  shade_p_value(obs_stat = Chisq_hat, direction = "greater")
+

Alternatively, visualizing the observed statistic using both of the null distributions,

-
visualize(null_distn, method = "both") +
-  shade_p_value(obs_stat = Chisq_hat, direction = "greater")
-

+
+visualize(null_dist, method = "both") +
+  shade_p_value(obs_stat = Chisq_hat, direction = "greater")
+

Note that the above code makes use of the randomization-based null distribution.

Calculating the p-value from the null distribution and observed statistic,

-
null_distn %>%
-  get_p_value(obs_stat = Chisq_hat, direction = "greater")
-
## # A tibble: 1 x 1
+
+null_dist %>%
+  get_p_value(obs_stat = Chisq_hat, direction = "greater")
+
## # A tibble: 1 × 1
 ##   p_value
 ##     <dbl>
-## 1   0.092
+## 1 0.108

Alternatively, using the wrapper to carry out the test,

-
gss %>%
-  chisq_test(formula = finrela ~ sex)
-
## # A tibble: 1 x 3
+
+gss %>%
+  chisq_test(formula = finrela ~ sex)
+
## # A tibble: 1 × 3
 ##   statistic chisq_df p_value
 ##       <dbl>    <int>   <dbl>
 ## 1      9.11        5   0.105

-One numerical variable, one categorical (2 levels) (diff in means)

+One numerical variable, one categorical (2 levels) (diff in means)

Calculating the observed statistic,

-
d_hat <- gss %>%
-  specify(age ~ college) %>%
-  calculate(stat = "diff in means", order = c("degree", "no degree"))
+
+d_hat <- gss %>% 
+  specify(age ~ college) %>% 
+  calculate(stat = "diff in means", order = c("degree", "no degree"))
+

Alternatively, using the observe() wrapper to calculate the observed statistic,

+
+d_hat <- gss %>% 
+  observe(age ~ college,
+          stat = "diff in means", order = c("degree", "no degree"))

Then, generating the null distribution,

-
null_distn <- gss %>%
-  specify(age ~ college) %>%
-  hypothesize(null = "independence") %>%
-  generate(reps = 1000, type = "permute") %>%
-  calculate(stat = "diff in means", order = c("degree", "no degree"))
+
+null_dist <- gss %>%
+  specify(age ~ college) %>%
+  hypothesize(null = "independence") %>%
+  generate(reps = 1000, type = "permute") %>%
+  calculate(stat = "diff in means", order = c("degree", "no degree"))

Visualizing the observed statistic alongside the null distribution,

-
visualize(null_distn) +
-  shade_p_value(obs_stat = d_hat, direction = "two-sided")
-

+
+visualize(null_dist) +
+  shade_p_value(obs_stat = d_hat, direction = "two-sided")
+

Calculating the p-value from the null distribution and observed statistic,

-
null_distn %>%
-  get_p_value(obs_stat = d_hat, direction = "two-sided")
-
## # A tibble: 1 x 1
+
+null_dist %>%
+  get_p_value(obs_stat = d_hat, direction = "two-sided")
+
## # A tibble: 1 × 1
 ##   p_value
 ##     <dbl>
-## 1   0.466
+## 1 0.438

-One numerical variable, one categorical (2 levels) (t)

+One numerical variable, one categorical (2 levels) (t)

Finding the standardized observed statistic,

-
t_hat <- gss %>%
-  specify(age ~ college) %>%
-  calculate(stat = "t", order = c("degree", "no degree"))
+
+t_hat <- gss %>% 
+  specify(age ~ college) %>% 
+  hypothesize(null = "independence") %>%
+  calculate(stat = "t", order = c("degree", "no degree"))
+

Alternatively, using the observe() wrapper to calculate the observed statistic,

+
+t_hat <- gss %>% 
+  observe(age ~ college,
+          stat = "t", order = c("degree", "no degree"))

Then, generating the null distribution,

-
null_distn <- gss %>%
-  specify(age ~ college) %>%
-  hypothesize(null = "independence") %>%
-  generate(reps = 1000, type = "permute") %>%
-  calculate(stat = "t", order = c("degree", "no degree"))
-

Alternatively, finding the null distribution using theoretical methods by skipping the generate() step,

-
null_distn_theoretical <- gss %>%
-  specify(age ~ college) %>%
-  hypothesize(null = "independence") %>%
-  calculate(stat = "t", order = c("degree", "no degree"))
+
+null_dist <- gss %>%
+  specify(age ~ college) %>%
+  hypothesize(null = "independence") %>%
+  generate(reps = 1000, type = "permute") %>%
+  calculate(stat = "t", order = c("degree", "no degree"))
+

Alternatively, finding the null distribution using theoretical methods using the assume() verb,

+
+null_dist_theory <- gss %>%
+  specify(age ~ college) %>%
+  assume("t")

Visualizing the observed statistic alongside the null distribution,

-
visualize(null_distn) +
-  shade_p_value(obs_stat = t_hat, direction = "two-sided")
-

+
+visualize(null_dist) +
+  shade_p_value(obs_stat = t_hat, direction = "two-sided")
+

Alternatively, visualizing the observed statistic using the theory-based null distribution,

-
visualize(null_distn_theoretical, method = "theoretical") +
-  shade_p_value(obs_stat = t_hat, direction = "two-sided")
-

+
+visualize(null_dist_theory) +
+  shade_p_value(obs_stat = t_hat, direction = "two-sided")
+

Alternatively, visualizing the observed statistic using both of the null distributions,

-
visualize(null_distn, method = "both") +
-  shade_p_value(obs_stat = t_hat, direction = "two-sided")
-

+
+visualize(null_dist, method = "both") +
+  shade_p_value(obs_stat = t_hat, direction = "two-sided")
+

Note that the above code makes use of the randomization-based null distribution.

Calculating the p-value from the null distribution and observed statistic,

-
null_distn %>%
-  get_p_value(obs_stat = t_hat, direction = "two-sided")
-
## # A tibble: 1 x 1
+
+null_dist %>%
+  get_p_value(obs_stat = t_hat, direction = "two-sided")
+
## # A tibble: 1 × 1
 ##   p_value
 ##     <dbl>
-## 1   0.408
+## 1 0.428

Note the similarities in this plot and the previous one.

-One numerical variable, one categorical (2 levels) (diff in medians)

+One numerical variable, one categorical (2 levels) (diff in medians)

Calculating the observed statistic,

-
d_hat <- gss %>%
-  specify(age ~ college) %>%
-  calculate(stat = "diff in medians", order = c("degree", "no degree"))
+
+d_hat <- gss %>% 
+  specify(age ~ college) %>% 
+  calculate(stat = "diff in medians", order = c("degree", "no degree"))
+

Alternatively, using the observe() wrapper to calculate the observed statistic,

+
+d_hat <- gss %>% 
+  observe(age ~ college,
+          stat = "diff in medians", order = c("degree", "no degree"))

Then, generating the null distribution,

-
null_distn <- gss %>%
-  specify(age ~ college) %>% # alt: response = age, explanatory = season
-  hypothesize(null = "independence") %>%
-  generate(reps = 1000, type = "permute") %>%
-  calculate(stat = "diff in medians", order = c("degree", "no degree"))
+
+null_dist <- gss %>%
+  specify(age ~ college) %>% # alt: response = age, explanatory = season
+  hypothesize(null = "independence") %>%
+  generate(reps = 1000, type = "permute") %>%
+  calculate(stat = "diff in medians", order = c("degree", "no degree"))

Visualizing the observed statistic alongside the null distribution,

-
visualize(null_distn) +
-  shade_p_value(obs_stat = d_hat, direction = "two-sided")
-

+
+visualize(null_dist) +
+  shade_p_value(obs_stat = d_hat, direction = "two-sided")
+

Calculating the p-value from the null distribution and observed statistic,

-
null_distn %>%
-  get_p_value(obs_stat = d_hat, direction = "two-sided")
-
## # A tibble: 1 x 1
+
+null_dist %>%
+  get_p_value(obs_stat = d_hat, direction = "two-sided")
+
## # A tibble: 1 × 1
 ##   p_value
 ##     <dbl>
-## 1   0.172
+## 1 0.194

-One numerical, one categorical (>2 levels) - ANOVA

+One numerical, one categorical (>2 levels) - ANOVA

Calculating the observed statistic,

-
F_hat <- gss %>%
-  specify(age ~ partyid) %>%
-  calculate(stat = "F")
+
+F_hat <- gss %>% 
+  specify(age ~ partyid) %>%
+  calculate(stat = "F")
+

Alternatively, using the observe() wrapper to calculate the observed statistic,

+
+F_hat <- gss %>% 
+  observe(age ~ partyid, stat = "F")

Then, generating the null distribution,

-
null_distn <- gss %>%
-   specify(age ~ partyid) %>%
-   hypothesize(null = "independence") %>%
-   generate(reps = 1000, type = "permute") %>%
-   calculate(stat = "F")
-

Alternatively, finding the null distribution using theoretical methods by skipping the generate() step,

-
null_distn_theoretical <- gss %>%
-   specify(age ~ partyid) %>%
-   hypothesize(null = "independence") %>%
-   calculate(stat = "F")
+
+null_dist <- gss %>%
+   specify(age ~ partyid) %>%
+   hypothesize(null = "independence") %>%
+   generate(reps = 1000, type = "permute") %>%
+   calculate(stat = "F")
+

Alternatively, finding the null distribution using theoretical methods using the assume() verb,

+
+null_dist_theory <- gss %>%
+   specify(age ~ partyid) %>%
+   hypothesize(null = "independence") %>%
+   assume(distribution = "F")

Visualizing the observed statistic alongside the null distribution,

-
visualize(null_distn) +
-  shade_p_value(obs_stat = F_hat, direction = "greater")
-

+
+visualize(null_dist) +
+  shade_p_value(obs_stat = F_hat, direction = "greater")
+

Alternatively, visualizing the observed statistic using the theory-based null distribution,

-
visualize(null_distn_theoretical, method = "theoretical") +
-  shade_p_value(obs_stat = F_hat, direction = "greater")
-

+
+visualize(null_dist_theory) +
+  shade_p_value(obs_stat = F_hat, direction = "greater")
+

Alternatively, visualizing the observed statistic using both of the null distributions,

-
visualize(null_distn, mdthod = "both") +
-  shade_p_value(obs_stat = F_hat, direction = "greater")
-

+
+visualize(null_dist, method = "both") +
+  shade_p_value(obs_stat = F_hat, direction = "greater")
+

Note that the above code makes use of the randomization-based null distribution.

Calculating the p-value from the null distribution and observed statistic,

-
null_distn %>%
-  get_p_value(obs_stat = F_hat, direction = "greater")
-
## # A tibble: 1 x 1
+
+null_dist %>%
+  get_p_value(obs_stat = F_hat, direction = "greater")
+
## # A tibble: 1 × 1
 ##   p_value
 ##     <dbl>
-## 1   0.046
+## 1 0.052

-Two numerical vars - SLR

+Two numerical vars - SLR

Calculating the observed statistic,

-
slope_hat <- gss %>%
-  specify(hours ~ age) %>%
-  calculate(stat = "slope")
+
+slope_hat <- gss %>% 
+  specify(hours ~ age) %>% 
+  calculate(stat = "slope")
+

Alternatively, using the observe() wrapper to calculate the observed statistic,

+
+slope_hat <- gss %>% 
+  observe(hours ~ age, stat = "slope")

Then, generating the null distribution,

-
null_distn <- gss %>%
-   specify(hours ~ age) %>%
-   hypothesize(null = "independence") %>%
-   generate(reps = 1000, type = "permute") %>%
-   calculate(stat = "slope")
+
+null_dist <- gss %>%
+   specify(hours ~ age) %>% 
+   hypothesize(null = "independence") %>%
+   generate(reps = 1000, type = "permute") %>%
+   calculate(stat = "slope")

Visualizing the observed statistic alongside the null distribution,

-
visualize(null_distn) +
-  shade_p_value(obs_stat = slope_hat, direction = "two-sided")
-

+
+visualize(null_dist) +
+  shade_p_value(obs_stat = slope_hat, direction = "two-sided")
+

Calculating the p-value from the null distribution and observed statistic,

-
null_distn %>%
-  get_p_value(obs_stat = slope_hat, direction = "two-sided")
-
## # A tibble: 1 x 1
+
+null_dist %>%
+  get_p_value(obs_stat = slope_hat, direction = "two-sided")
+
## # A tibble: 1 × 1
 ##   p_value
 ##     <dbl>
-## 1   0.838
+## 1 0.88

-Two numerical vars - correlation

+Two numerical vars - correlation

Calculating the observed statistic,

-
correlation_hat <- gss %>%
-  specify(hours ~ age) %>%
-  calculate(stat = "correlation")
+
+correlation_hat <- gss %>% 
+  specify(hours ~ age) %>% 
+  calculate(stat = "correlation")
+

Alternatively, using the observe() wrapper to calculate the observed statistic,

+
+correlation_hat <- gss %>% 
+  observe(hours ~ age, stat = "correlation")

Then, generating the null distribution,

-
null_distn <- gss %>%
-   specify(hours ~ age) %>%
-   hypothesize(null = "independence") %>%
-   generate(reps = 1000, type = "permute") %>%
-   calculate(stat = "correlation")
+
+null_dist <- gss %>%
+   specify(hours ~ age) %>% 
+   hypothesize(null = "independence") %>%
+   generate(reps = 1000, type = "permute") %>%
+   calculate(stat = "correlation")

Visualizing the observed statistic alongside the null distribution,

-
visualize(null_distn) +
-  shade_p_value(obs_stat = correlation_hat, direction = "two-sided")
-

+
+visualize(null_dist) +
+  shade_p_value(obs_stat = correlation_hat, direction = "two-sided")
+

Calculating the p-value from the null distribution and observed statistic,

-
null_distn %>%
-  get_p_value(obs_stat = correlation_hat, direction = "two-sided")
-
## # A tibble: 1 x 1
+
+null_dist %>%
+  get_p_value(obs_stat = correlation_hat, direction = "two-sided")
+
## # A tibble: 1 × 1
 ##   p_value
 ##     <dbl>
-## 1     0.9
+## 1 0.876

-Two numerical vars - SLR (t)

+Two numerical vars - SLR (t)

Not currently implemented since \(t\) could refer to standardized slope or standardized correlation.

+
+

+Multiple explanatory variables

+

Calculating the observed fit,

+
+obs_fit <- gss %>%
+  specify(hours ~ age + college) %>%
+  fit()
+

Generating a distribution of fits with the response variable permuted,

+
+null_dist <- gss %>%
+  specify(hours ~ age + college) %>%
+  hypothesize(null = "independence") %>%
+  generate(reps = 1000, type = "permute") %>%
+  fit()
+

Generating a distribution of fits where each explanatory variable is permuted independently,

+
+null_dist2 <- gss %>%
+  specify(hours ~ age + college) %>%
+  hypothesize(null = "independence") %>%
+  generate(reps = 1000, type = "permute", variables = c(age, college)) %>%
+  fit()
+

Visualizing the observed fit alongside the null fits,

+
+visualize(null_dist) +
+  shade_p_value(obs_stat = obs_fit, direction = "two-sided")
+

+

Calculating p-values from the null distribution and observed fit,

+
+null_dist %>%
+  get_p_value(obs_stat = obs_fit, direction = "two-sided")
+
## # A tibble: 3 × 2
+##   term          p_value
+##   <chr>           <dbl>
+## 1 age             0.852
+## 2 collegedegree   0.278
+## 3 intercept       0.634
+

Note that this fit()-based workflow can be applied to use cases with differing numbers of explanatory variables and explanatory variable types.

+

-Confidence intervals

+Confidence intervals

-One numerical (one mean)

+One numerical (one mean)

Finding the observed statistic,

-
x_bar <- gss %>%
-  specify(response = hours) %>%
-  calculate(stat = "mean")
-

Then, generating the null distribution,

-
boot <- gss %>%
-   specify(response = hours) %>%
-   generate(reps = 1000, type = "bootstrap") %>%
-   calculate(stat = "mean")
-

Use the null distribution to find a confidence interval,

-
percentile_ci <- get_ci(boot)
-

Visualizing the observed statistic alongside the null distribution,

-
visualize(boot) +
-  shade_confidence_interval(endpoints = percentile_ci)
-

-

Alternatively, use the null distribution to find a confidence interval using the standard error,

-
standard_error_ci <- get_ci(boot, type = "se", point_estimate = x_bar)
+
+x_bar <- gss %>% 
+  specify(response = hours) %>%
+  calculate(stat = "mean")
+

Alternatively, using the observe() wrapper to calculate the observed statistic,

+
+x_bar <- gss %>% 
+  observe(response = hours, stat = "mean")
+

Then, generating a bootstrap distribution,

+
+boot_dist <- gss %>%
+   specify(response = hours) %>%
+   generate(reps = 1000, type = "bootstrap") %>%
+   calculate(stat = "mean")
+

Use the bootstrap distribution to find a confidence interval,

+
+percentile_ci <- get_ci(boot_dist)
+

Visualizing the observed statistic alongside the distribution,

+
+visualize(boot_dist) +
+  shade_confidence_interval(endpoints = percentile_ci)
+

+

Alternatively, use the bootstrap distribution to find a confidence interval using the standard error,

+
+standard_error_ci <- get_ci(boot_dist, type = "se", point_estimate = x_bar)
 
-visualize(boot) +
-  shade_confidence_interval(endpoints = standard_error_ci)
-

+visualize(boot_dist) + + shade_confidence_interval(endpoints = standard_error_ci)
+

+

Instead of a simulation-based bootstrap distribution, we can also define a theory-based sampling distribution,

+
+sampling_dist <- gss %>%
+   specify(response = hours) %>%
+   assume(distribution = "t")
+

Visualization and calculation of confidence intervals interfaces in the same way as with the simulation-based distribution,

+
+theor_ci <- get_ci(sampling_dist, point_estimate = x_bar)
+
+theor_ci
+
## # A tibble: 1 × 2
+##   lower_ci upper_ci
+##      <dbl>    <dbl>
+## 1     40.1     42.7
+
+visualize(sampling_dist) +
+  shade_confidence_interval(endpoints = theor_ci)
+

+

Note that the t distribution is recentered and rescaled to lie on the scale of the observed data. infer does not support confidence intervals on means via the z distribution.

-One numerical (one mean - standardized)

+One numerical (one mean - standardized)

Finding the observed statistic,

-
t_hat <- gss %>%
-  specify(response = hours) %>%
-  hypothesize(null = "point", mu = 40) %>%
-  calculate(stat = "t")
-

Then, generating the null distribution,

-
boot <- gss %>%
-   specify(response = hours) %>%
-   hypothesize(null = "point", mu = 40) %>%
-   generate(reps = 1000, type = "bootstrap") %>%
-   calculate(stat = "t")
-

Use the null distribution to find a confidence interval,

-
percentile_ci <- get_ci(boot)
-

Visualizing the observed statistic alongside the null distribution,

-
visualize(boot) +
-  shade_confidence_interval(endpoints = percentile_ci)
-

-

Alternatively, use the null distribution to find a confidence interval using the standard error,

-
standard_error_ci <- boot %>%
-  get_ci(type = "se", point_estimate = t_hat)
+
+t_hat <- gss %>% 
+  specify(response = hours) %>%
+  hypothesize(null = "point", mu = 40) %>%
+  calculate(stat = "t")
+

Alternatively, using the observe() wrapper to calculate the observed statistic,

+
+t_hat <- gss %>% 
+  observe(response = hours,
+          null = "point", mu = 40,
+          stat = "t")
+

Then, generating the bootstrap distribution,

+
+boot_dist <- gss %>%
+   specify(response = hours) %>%
+   generate(reps = 1000, type = "bootstrap") %>%
+   calculate(stat = "t")
+

Use the bootstrap distribution to find a confidence interval,

+
+percentile_ci <- get_ci(boot_dist)
+

Visualizing the observed statistic alongside the distribution,

+
+visualize(boot_dist) +
+  shade_confidence_interval(endpoints = percentile_ci)
+

+

Alternatively, use the bootstrap distribution to find a confidence interval using the standard error,

+
+standard_error_ci <- boot_dist %>%
+  get_ci(type = "se", point_estimate = t_hat)
 
-visualize(boot) +
-  shade_confidence_interval(endpoints = standard_error_ci)
-

+visualize(boot_dist) + + shade_confidence_interval(endpoints = standard_error_ci)
+

+

See the above subsection (one mean) for a theory-based approach. Note that infer does not support confidence intervals on means via the z distribution.

-One categorical (one proportion)

+One categorical (one proportion)

Finding the observed statistic,

-
p_hat <- gss %>%
-   specify(response = sex, success = "female") %>%
-   calculate(stat = "prop")
-

Then, generating the null distribution,

-
boot <- gss %>%
- specify(response = sex, success = "female") %>%
- generate(reps = 1000, type = "bootstrap") %>%
- calculate(stat = "prop")
-

Use the null distribution to find a confidence interval,

-
percentile_ci <- get_ci(boot)
-

Visualizing the observed statistic alongside the null distribution,

-
visualize(boot) +
-  shade_confidence_interval(endpoints = percentile_ci)
-

-

Alternatively, use the null distribution to find a confidence interval using the standard error,

-
standard_error_ci <- boot %>%
-  get_ci(type = "se", point_estimate = p_hat)
+
+p_hat <- gss %>% 
+   specify(response = sex, success = "female") %>%
+   calculate(stat = "prop")
+

Alternatively, using the observe() wrapper to calculate the observed statistic,

+
+p_hat <- gss %>% 
+   observe(response = sex, success = "female", stat = "prop")
+

Then, generating a bootstrap distribution,

+
+boot_dist <- gss %>%
+ specify(response = sex, success = "female") %>%
+ generate(reps = 1000, type = "bootstrap") %>%
+ calculate(stat = "prop")
+

Use the bootstrap distribution to find a confidence interval,

+
+percentile_ci <- get_ci(boot_dist)
+

Visualizing the observed statistic alongside the distribution,

+
+visualize(boot_dist) +
+  shade_confidence_interval(endpoints = percentile_ci)
+

+

Alternatively, use the bootstrap distribution to find a confidence interval using the standard error,

+
+standard_error_ci <- boot_dist %>%
+  get_ci(type = "se", point_estimate = p_hat)
+
+visualize(boot_dist) +
+  shade_confidence_interval(endpoints = standard_error_ci)
+

+

Instead of a simulation-based bootstrap distribution, we can also define a theory-based sampling distribution,

+
+sampling_dist <- gss %>%
+   specify(response = sex, success = "female") %>%
+   assume(distribution = "z")
+

Visualization and calculation of confidence intervals interfaces in the same way as with the simulation-based distribution,

+
+theor_ci <- get_ci(sampling_dist, point_estimate = p_hat)
 
-visualize(boot) +
-  shade_confidence_interval(endpoints = standard_error_ci)
-

+theor_ci
+
## # A tibble: 1 × 2
+##   lower_ci upper_ci
+##      <dbl>    <dbl>
+## 1    0.430    0.518
+
+visualize(sampling_dist) +
+  shade_confidence_interval(endpoints = theor_ci)
+

+

Note that the z distribution is recentered and rescaled to lie on the scale of the observed data. infer does not support confidence intervals on means via the z distribution.

-One categorical variable (standardized proportion \(z\))

-

Not yet implemented.

+One categorical variable (standardized proportion \(z\)) +

See the above subsection (one proportion) for a theory-based approach.

-One numerical variable, one categorical (2 levels) (diff in means)

+One numerical variable, one categorical (2 levels) (diff in means)

Finding the observed statistic,

-
d_hat <- gss %>%
-  specify(hours ~ college) %>%
-  calculate(stat = "diff in means", order = c("degree", "no degree"))
-

Then, generating the null distribution,

-
boot <- gss %>%
-   specify(hours ~ college) %>%
-   generate(reps = 1000, type = "bootstrap") %>%
-   calculate(stat = "diff in means", order = c("degree", "no degree"))
-

Use the null distribution to find a confidence interval,

-
percentile_ci <- get_ci(boot)
-

Visualizing the observed statistic alongside the null distribution,

-
visualize(boot) +
-  shade_confidence_interval(endpoints = percentile_ci)
-

-

Alternatively, use the null distribution to find a confidence interval using the standard error,

-
standard_error_ci <- boot %>%
-  get_ci(type = "se", point_estimate = d_hat)
+
+d_hat <- gss %>%
+  specify(hours ~ college) %>%
+  calculate(stat = "diff in means", order = c("degree", "no degree"))
+

Alternatively, using the observe() wrapper to calculate the observed statistic,

+
+d_hat <- gss %>%
+  observe(hours ~ college,
+          stat = "diff in means", order = c("degree", "no degree"))
+

Then, generating a bootstrap distribution,

+
+boot_dist <- gss %>%
+   specify(hours ~ college) %>%
+   generate(reps = 1000, type = "bootstrap") %>%
+   calculate(stat = "diff in means", order = c("degree", "no degree"))
+

Use the bootstrap distribution to find a confidence interval,

+
+percentile_ci <- get_ci(boot_dist)
+

Visualizing the observed statistic alongside the distribution,

+
+visualize(boot_dist) +
+  shade_confidence_interval(endpoints = percentile_ci)
+

+

Alternatively, use the bootstrap distribution to find a confidence interval using the standard error,

+
+standard_error_ci <- boot_dist %>%
+  get_ci(type = "se", point_estimate = d_hat)
 
-visualize(boot) +
-  shade_confidence_interval(endpoints = standard_error_ci)
-

+visualize(boot_dist) + + shade_confidence_interval(endpoints = standard_error_ci)
+

+

Instead of a simulation-based bootstrap distribution, we can also define a theory-based sampling distribution,

+
+sampling_dist <- gss %>%
+   specify(hours ~ college) %>%
+   assume(distribution = "t")
+

Visualization and calculation of confidence intervals interfaces in the same way as with the simulation-based distribution,

+
+theor_ci <- get_ci(sampling_dist, point_estimate = d_hat)
+
+theor_ci
+
## # A tibble: 1 × 2
+##   lower_ci upper_ci
+##      <dbl>    <dbl>
+## 1    -1.16     4.24
+
+visualize(sampling_dist) +
+  shade_confidence_interval(endpoints = theor_ci)
+

+

Note that the t distribution is recentered and rescaled to lie on the scale of the observed data.

-One numerical variable, one categorical (2 levels) (t)

+One numerical variable, one categorical (2 levels) (t)

Finding the standardized point estimate,

-
t_hat <- gss %>%
-  specify(hours ~ college) %>%
-  calculate(stat = "t", order = c("degree", "no degree"))
-

Then, generating the null distribution,

-
boot <- gss %>%
-   specify(hours ~ college) %>%
-   generate(reps = 1000, type = "bootstrap") %>%
-   calculate(stat = "t", order = c("degree", "no degree"))
-

Use the null distribution to find a confidence interval,

-
percentile_ci <- get_ci(boot)
-

Visualizing the observed statistic alongside the null distribution,

-
visualize(boot) +
-  shade_confidence_interval(endpoints = percentile_ci)
-

-

Alternatively, use the null distribution to find a confidence interval using the standard error,

-
standard_error_ci <- boot %>%
-  get_ci(type = "se", point_estimate = t_hat)
+
+t_hat <- gss %>%
+  specify(hours ~ college) %>%
+  calculate(stat = "t", order = c("degree", "no degree"))
+

Alternatively, using the observe() wrapper to calculate the observed statistic,

+
+t_hat <- gss %>%
+  observe(hours ~ college,
+          stat = "t", order = c("degree", "no degree"))
+

Then, generating a bootstrap distribution,

+
+boot_dist <- gss %>%
+   specify(hours ~ college) %>%
+   generate(reps = 1000, type = "bootstrap") %>%
+   calculate(stat = "t", order = c("degree", "no degree"))
+

Use the bootstrap distribution to find a confidence interval,

+
+percentile_ci <- get_ci(boot_dist)
+

Visualizing the observed statistic alongside the distribution,

+
+visualize(boot_dist) +
+  shade_confidence_interval(endpoints = percentile_ci)
+

+

Alternatively, use the bootstrap distribution to find a confidence interval using the standard error,

+
+standard_error_ci <- boot_dist %>%
+  get_ci(type = "se", point_estimate = t_hat)
 
-visualize(boot) +
-  shade_confidence_interval(endpoints = standard_error_ci)
-

+visualize(boot_dist) + + shade_confidence_interval(endpoints = standard_error_ci)
+

+

See the above subsection (diff in means) for a theory-based approach. infer does not support confidence intervals on means via the z distribution.

-Two categorical variables (diff in proportions)

+Two categorical variables (diff in proportions)

Finding the observed statistic,

-
d_hat <- gss %>%
-  specify(college ~ sex, success = "degree") %>%
-  calculate(stat = "diff in props", order = c("female", "male"))
-

Then, generating the null distribution,

-
boot <- gss %>%
-  specify(college ~ sex, success = "degree") %>%
-  generate(reps = 1000, type = "bootstrap") %>%
-  calculate(stat = "diff in props", order = c("female", "male"))
-

Use the null distribution to find a confidence interval,

-
percentile_ci <- get_ci(boot)
-

Visualizing the observed statistic alongside the null distribution,

-
visualize(boot) +
-  shade_confidence_interval(endpoints = percentile_ci)
-

-

Alternatively, use the null distribution to find a confidence interval using the standard error,

-
standard_error_ci <- boot %>%
-  get_ci(type = "se", point_estimate = d_hat)
-
-visualize(boot) +
-  shade_confidence_interval(endpoints = standard_error_ci)
-

+
+d_hat <- gss %>% 
+  specify(college ~ sex, success = "degree") %>%
+  calculate(stat = "diff in props", order = c("female", "male"))
+

Alternatively, using the observe() wrapper to calculate the observed statistic,

+
+d_hat <- gss %>% 
+  observe(college ~ sex, success = "degree",
+          stat = "diff in props", order = c("female", "male"))
+

Then, generating a bootstrap distribution,

+
+boot_dist <- gss %>%
+  specify(college ~ sex, success = "degree") %>%
+  generate(reps = 1000, type = "bootstrap") %>% 
+  calculate(stat = "diff in props", order = c("female", "male"))
+

Use the bootstrap distribution to find a confidence interval,

+
+percentile_ci <- get_ci(boot_dist)
+

Visualizing the observed statistic alongside the distribution,

+
+visualize(boot_dist) +
+  shade_confidence_interval(endpoints = percentile_ci)
+

+

Alternatively, use the bootstrap distribution to find a confidence interval using the standard error,

+
+standard_error_ci <- boot_dist %>%
+  get_ci(type = "se", point_estimate = d_hat)
+
+visualize(boot_dist) +
+  shade_confidence_interval(endpoints = standard_error_ci)
+

+

Instead of a simulation-based bootstrap distribution, we can also define a theory-based sampling distribution,

+
+sampling_dist <- gss %>% 
+  specify(college ~ sex, success = "degree") %>%
+   assume(distribution = "z")
+

Visualization and calculation of confidence intervals interfaces in the same way as with the simulation-based distribution,

+
+theor_ci <- get_ci(sampling_dist, point_estimate = d_hat)
+
+theor_ci
+
## # A tibble: 1 × 2
+##   lower_ci upper_ci
+##      <dbl>    <dbl>
+## 1  -0.0794   0.0878
+
+visualize(sampling_dist) +
+  shade_confidence_interval(endpoints = theor_ci)
+

+

Note that the z distribution is recentered and rescaled to lie on the scale of the observed data.

-Two categorical variables (z)

+Two categorical variables (z)

Finding the standardized point estimate,

-
z_hat <- gss %>%
-  specify(college ~ sex, success = "degree") %>%
-  calculate(stat = "z", order = c("female", "male"))
-

Then, generating the null distribution,

-
boot <- gss %>%
-  specify(college ~ sex, success = "degree") %>%
-  generate(reps = 1000, type = "bootstrap") %>%
-  calculate(stat = "z", order = c("female", "male"))
-

Use the null distribution to find a confidence interval,

-
percentile_ci <- get_ci(boot)
-

Visualizing the observed statistic alongside the null distribution,

-
visualize(boot) +
-  shade_confidence_interval(endpoints = percentile_ci)
-

-

Alternatively, use the null distribution to find a confidence interval using the standard error,

-
standard_error_ci <- boot %>%
-  get_ci(type = "se", point_estimate = z_hat)
-
-visualize(boot) +
-  shade_confidence_interval(endpoints = standard_error_ci)
-

+
+z_hat <- gss %>% 
+  specify(college ~ sex, success = "degree") %>%
+  calculate(stat = "z", order = c("female", "male"))
+

Alternatively, using the observe() wrapper to calculate the observed statistic,

+
+z_hat <- gss %>% 
+  observe(college ~ sex, success = "degree",
+          stat = "z", order = c("female", "male"))
+

Then, generating a bootstrap distribution,

+
+boot_dist <- gss %>%
+  specify(college ~ sex, success = "degree") %>%
+  generate(reps = 1000, type = "bootstrap") %>% 
+  calculate(stat = "z", order = c("female", "male"))
+

Use the bootstrap distribution to find a confidence interval,

+
+percentile_ci <- get_ci(boot_dist)
+

Visualizing the observed statistic alongside the distribution,

+
+visualize(boot_dist) +
+  shade_confidence_interval(endpoints = percentile_ci)
+

+

Alternatively, use the bootstrap distribution to find a confidence interval using the standard error,

+
+standard_error_ci <- boot_dist %>%
+  get_ci(type = "se", point_estimate = z_hat)
+
+visualize(boot_dist) +
+  shade_confidence_interval(endpoints = standard_error_ci)
+

+

See the above subsection (diff in props) for a theory-based approach.

-Two numerical vars - SLR

+Two numerical vars - SLR

Finding the observed statistic,

-
slope_hat <- gss %>%
-  specify(hours ~ age) %>%
-  calculate(stat = "slope")
-

Then, generating the null distribution,

-
boot <- gss %>%
-   specify(hours ~ age) %>%
-   generate(reps = 1000, type = "bootstrap") %>%
-   calculate(stat = "slope")
-

Use the null distribution to find a confidence interval,

-
percentile_ci <- get_ci(boot)
-

Visualizing the observed statistic alongside the null distribution,

-
visualize(boot) +
-  shade_confidence_interval(endpoints = percentile_ci)
-

-

Alternatively, use the null distribution to find a confidence interval using the standard error,

-
standard_error_ci <- boot %>%
-  get_ci(type = "se", point_estimate = slope_hat)
+
+slope_hat <- gss %>% 
+  specify(hours ~ age) %>%
+  calculate(stat = "slope")
+

Alternatively, using the observe() wrapper to calculate the observed statistic,

+
+slope_hat <- gss %>% 
+  observe(hours ~ age, stat = "slope")
+

Then, generating a bootstrap distribution,

+
+boot_dist <- gss %>%
+   specify(hours ~ age) %>% 
+   generate(reps = 1000, type = "bootstrap") %>%
+   calculate(stat = "slope")
+

Use the bootstrap distribution to find a confidence interval,

+
+percentile_ci <- get_ci(boot_dist)
+

Visualizing the observed statistic alongside the distribution,

+
+visualize(boot_dist) +
+  shade_confidence_interval(endpoints = percentile_ci)
+

+

Alternatively, use the bootstrap distribution to find a confidence interval using the standard error,

+
+standard_error_ci <- boot_dist %>%
+  get_ci(type = "se", point_estimate = slope_hat)
 
-visualize(boot) +
-  shade_confidence_interval(endpoints = standard_error_ci)
-

+visualize(boot_dist) + + shade_confidence_interval(endpoints = standard_error_ci)
+

-Two numerical vars - correlation

+Two numerical vars - correlation

Finding the observed statistic,

-
correlation_hat <- gss %>%
-  specify(hours ~ age) %>%
-  calculate(stat = "correlation")
-

Then, generating the null distribution,

-
boot <- gss %>%
-   specify(hours ~ age) %>%
-   generate(reps = 1000, type = "bootstrap") %>%
-   calculate(stat = "correlation")
-

Use the null distribution to find a confidence interval,

-
percentile_ci <- get_ci(boot)
-

Visualizing the observed statistic alongside the null distribution,

-
visualize(boot) +
-  shade_confidence_interval(endpoints = percentile_ci)
-

-

Alternatively, use the null distribution to find a confidence interval using the standard error,

-
standard_error_ci <- boot %>%
-  get_ci(type = "se", point_estimate = correlation_hat)
+
+correlation_hat <- gss %>% 
+  specify(hours ~ age) %>%
+  calculate(stat = "correlation")
+

Alternatively, using the observe() wrapper to calculate the observed statistic,

+
+correlation_hat <- gss %>% 
+  observe(hours ~ age, stat = "correlation")
+

Then, generating a bootstrap distribution,

+
+boot_dist <- gss %>%
+   specify(hours ~ age) %>% 
+   generate(reps = 1000, type = "bootstrap") %>%
+   calculate(stat = "correlation")
+

Use the bootstrap distribution to find a confidence interval,

+
+percentile_ci <- get_ci(boot_dist)
+

Visualizing the observed statistic alongside the distribution,

+
+visualize(boot_dist) +
+  shade_confidence_interval(endpoints = percentile_ci)
+

+

Alternatively, use the bootstrap distribution to find a confidence interval using the standard error,

+
+standard_error_ci <- boot_dist %>%
+  get_ci(type = "se", point_estimate = correlation_hat)
 
-visualize(boot) +
-  shade_confidence_interval(endpoints = standard_error_ci)
-

+visualize(boot_dist) + + shade_confidence_interval(endpoints = standard_error_ci)
+

-Two numerical vars - t

+Two numerical vars - t

Not currently implemented since \(t\) could refer to standardized slope or standardized correlation.

+
+

+Multiple explanatory variables

+

Calculating the observed fit,

+
+obs_fit <- gss %>%
+  specify(hours ~ age + college) %>%
+  fit()
+

Generating a distribution of fits with the response variable permuted,

+
+null_dist <- gss %>%
+  specify(hours ~ age + college) %>%
+  hypothesize(null = "independence") %>%
+  generate(reps = 1000, type = "permute") %>%
+  fit()
+

Alternatively, generating a distribution of fits where each explanatory variable is permuted independently,

+
+null_dist2 <- gss %>%
+  specify(hours ~ age + college) %>%
+  hypothesize(null = "independence") %>%
+  generate(reps = 1000, type = "permute", variables = c(age, college)) %>%
+  fit()
+

Calculating confidence intervals from the null fits,

+
+conf_ints <- 
+  get_confidence_interval(
+    null_dist, 
+    level = .95, 
+    point_estimate = obs_fit
+  )
+

Visualizing the observed fit alongside the null fits,

+
+visualize(null_dist) +
+  shade_confidence_interval(endpoints = conf_ints)
+

+

Note that this fit()-based workflow can be applied to use cases with differing numbers of explanatory variables and explanatory variable types.

+
@@ -983,32 +1400,24 @@

-
- - + + + + + diff --git a/docs/articles/observed_stat_examples_files/figure-html/unnamed-chunk-101-1.png b/docs/articles/observed_stat_examples_files/figure-html/unnamed-chunk-101-1.png new file mode 100644 index 00000000..6d8c2ab3 Binary files /dev/null and b/docs/articles/observed_stat_examples_files/figure-html/unnamed-chunk-101-1.png differ diff --git a/docs/articles/observed_stat_examples_files/figure-html/unnamed-chunk-106-1.png b/docs/articles/observed_stat_examples_files/figure-html/unnamed-chunk-106-1.png new file mode 100644 index 00000000..8cffd218 Binary files /dev/null and b/docs/articles/observed_stat_examples_files/figure-html/unnamed-chunk-106-1.png differ diff --git a/docs/articles/observed_stat_examples_files/figure-html/unnamed-chunk-11-1.png b/docs/articles/observed_stat_examples_files/figure-html/unnamed-chunk-11-1.png index 5be01f7d..dda62af4 100644 Binary files a/docs/articles/observed_stat_examples_files/figure-html/unnamed-chunk-11-1.png and b/docs/articles/observed_stat_examples_files/figure-html/unnamed-chunk-11-1.png differ diff --git a/docs/articles/observed_stat_examples_files/figure-html/unnamed-chunk-111-1.png b/docs/articles/observed_stat_examples_files/figure-html/unnamed-chunk-111-1.png new file mode 100644 index 00000000..219442ca Binary files /dev/null and b/docs/articles/observed_stat_examples_files/figure-html/unnamed-chunk-111-1.png differ diff --git a/docs/articles/observed_stat_examples_files/figure-html/unnamed-chunk-117-1.png b/docs/articles/observed_stat_examples_files/figure-html/unnamed-chunk-117-1.png new file mode 100644 index 00000000..8fb3dda1 Binary files /dev/null and b/docs/articles/observed_stat_examples_files/figure-html/unnamed-chunk-117-1.png differ diff --git a/docs/articles/observed_stat_examples_files/figure-html/unnamed-chunk-118-1.png b/docs/articles/observed_stat_examples_files/figure-html/unnamed-chunk-118-1.png new file mode 100644 index 00000000..455f6d14 Binary files /dev/null and b/docs/articles/observed_stat_examples_files/figure-html/unnamed-chunk-118-1.png differ diff --git a/docs/articles/observed_stat_examples_files/figure-html/unnamed-chunk-12-1.png b/docs/articles/observed_stat_examples_files/figure-html/unnamed-chunk-12-1.png index 31a2b924..a096db3b 100644 Binary files a/docs/articles/observed_stat_examples_files/figure-html/unnamed-chunk-12-1.png and b/docs/articles/observed_stat_examples_files/figure-html/unnamed-chunk-12-1.png differ diff --git a/docs/articles/observed_stat_examples_files/figure-html/unnamed-chunk-120-1.png b/docs/articles/observed_stat_examples_files/figure-html/unnamed-chunk-120-1.png index d31770d2..98c6c7cd 100644 Binary files a/docs/articles/observed_stat_examples_files/figure-html/unnamed-chunk-120-1.png and b/docs/articles/observed_stat_examples_files/figure-html/unnamed-chunk-120-1.png differ diff --git a/docs/articles/observed_stat_examples_files/figure-html/unnamed-chunk-125-1.png b/docs/articles/observed_stat_examples_files/figure-html/unnamed-chunk-125-1.png index ece8eea0..848fba82 100644 Binary files a/docs/articles/observed_stat_examples_files/figure-html/unnamed-chunk-125-1.png and b/docs/articles/observed_stat_examples_files/figure-html/unnamed-chunk-125-1.png differ diff --git a/docs/articles/observed_stat_examples_files/figure-html/unnamed-chunk-126-1.png b/docs/articles/observed_stat_examples_files/figure-html/unnamed-chunk-126-1.png new file mode 100644 index 00000000..83e76a9a Binary files /dev/null and b/docs/articles/observed_stat_examples_files/figure-html/unnamed-chunk-126-1.png differ diff --git a/docs/articles/observed_stat_examples_files/figure-html/unnamed-chunk-13-1.png b/docs/articles/observed_stat_examples_files/figure-html/unnamed-chunk-13-1.png index f2dab0df..9d7489cb 100644 Binary files a/docs/articles/observed_stat_examples_files/figure-html/unnamed-chunk-13-1.png and b/docs/articles/observed_stat_examples_files/figure-html/unnamed-chunk-13-1.png differ diff --git a/docs/articles/observed_stat_examples_files/figure-html/unnamed-chunk-131-1.png b/docs/articles/observed_stat_examples_files/figure-html/unnamed-chunk-131-1.png new file mode 100644 index 00000000..65a31147 Binary files /dev/null and b/docs/articles/observed_stat_examples_files/figure-html/unnamed-chunk-131-1.png differ diff --git a/docs/articles/observed_stat_examples_files/figure-html/unnamed-chunk-132-1.png b/docs/articles/observed_stat_examples_files/figure-html/unnamed-chunk-132-1.png new file mode 100644 index 00000000..61b838b4 Binary files /dev/null and b/docs/articles/observed_stat_examples_files/figure-html/unnamed-chunk-132-1.png differ diff --git a/docs/articles/observed_stat_examples_files/figure-html/unnamed-chunk-134-1.png b/docs/articles/observed_stat_examples_files/figure-html/unnamed-chunk-134-1.png index 93fdeeba..d389c213 100644 Binary files a/docs/articles/observed_stat_examples_files/figure-html/unnamed-chunk-134-1.png and b/docs/articles/observed_stat_examples_files/figure-html/unnamed-chunk-134-1.png differ diff --git a/docs/articles/observed_stat_examples_files/figure-html/unnamed-chunk-139-1.png b/docs/articles/observed_stat_examples_files/figure-html/unnamed-chunk-139-1.png new file mode 100644 index 00000000..3a847aa8 Binary files /dev/null and b/docs/articles/observed_stat_examples_files/figure-html/unnamed-chunk-139-1.png differ diff --git a/docs/articles/observed_stat_examples_files/figure-html/unnamed-chunk-140-1.png b/docs/articles/observed_stat_examples_files/figure-html/unnamed-chunk-140-1.png new file mode 100644 index 00000000..535e07ef Binary files /dev/null and b/docs/articles/observed_stat_examples_files/figure-html/unnamed-chunk-140-1.png differ diff --git a/docs/articles/observed_stat_examples_files/figure-html/unnamed-chunk-142-1.png b/docs/articles/observed_stat_examples_files/figure-html/unnamed-chunk-142-1.png new file mode 100644 index 00000000..265c2347 Binary files /dev/null and b/docs/articles/observed_stat_examples_files/figure-html/unnamed-chunk-142-1.png differ diff --git a/docs/articles/observed_stat_examples_files/figure-html/unnamed-chunk-147-1.png b/docs/articles/observed_stat_examples_files/figure-html/unnamed-chunk-147-1.png new file mode 100644 index 00000000..779386f0 Binary files /dev/null and b/docs/articles/observed_stat_examples_files/figure-html/unnamed-chunk-147-1.png differ diff --git a/docs/articles/observed_stat_examples_files/figure-html/unnamed-chunk-148-1.png b/docs/articles/observed_stat_examples_files/figure-html/unnamed-chunk-148-1.png new file mode 100644 index 00000000..6ddef50b Binary files /dev/null and b/docs/articles/observed_stat_examples_files/figure-html/unnamed-chunk-148-1.png differ diff --git a/docs/articles/observed_stat_examples_files/figure-html/unnamed-chunk-153-1.png b/docs/articles/observed_stat_examples_files/figure-html/unnamed-chunk-153-1.png new file mode 100644 index 00000000..dbaa1024 Binary files /dev/null and b/docs/articles/observed_stat_examples_files/figure-html/unnamed-chunk-153-1.png differ diff --git a/docs/articles/observed_stat_examples_files/figure-html/unnamed-chunk-154-1.png b/docs/articles/observed_stat_examples_files/figure-html/unnamed-chunk-154-1.png new file mode 100644 index 00000000..64cbd9fb Binary files /dev/null and b/docs/articles/observed_stat_examples_files/figure-html/unnamed-chunk-154-1.png differ diff --git a/docs/articles/observed_stat_examples_files/figure-html/unnamed-chunk-156-1.png b/docs/articles/observed_stat_examples_files/figure-html/unnamed-chunk-156-1.png new file mode 100644 index 00000000..6cf3dde2 Binary files /dev/null and b/docs/articles/observed_stat_examples_files/figure-html/unnamed-chunk-156-1.png differ diff --git a/docs/articles/observed_stat_examples_files/figure-html/unnamed-chunk-161-1.png b/docs/articles/observed_stat_examples_files/figure-html/unnamed-chunk-161-1.png new file mode 100644 index 00000000..1caa66ea Binary files /dev/null and b/docs/articles/observed_stat_examples_files/figure-html/unnamed-chunk-161-1.png differ diff --git a/docs/articles/observed_stat_examples_files/figure-html/unnamed-chunk-162-1.png b/docs/articles/observed_stat_examples_files/figure-html/unnamed-chunk-162-1.png new file mode 100644 index 00000000..3c179f7f Binary files /dev/null and b/docs/articles/observed_stat_examples_files/figure-html/unnamed-chunk-162-1.png differ diff --git a/docs/articles/observed_stat_examples_files/figure-html/unnamed-chunk-167-1.png b/docs/articles/observed_stat_examples_files/figure-html/unnamed-chunk-167-1.png new file mode 100644 index 00000000..82fdc1da Binary files /dev/null and b/docs/articles/observed_stat_examples_files/figure-html/unnamed-chunk-167-1.png differ diff --git a/docs/articles/observed_stat_examples_files/figure-html/unnamed-chunk-168-1.png b/docs/articles/observed_stat_examples_files/figure-html/unnamed-chunk-168-1.png new file mode 100644 index 00000000..05220ce4 Binary files /dev/null and b/docs/articles/observed_stat_examples_files/figure-html/unnamed-chunk-168-1.png differ diff --git a/docs/articles/observed_stat_examples_files/figure-html/unnamed-chunk-173-1.png b/docs/articles/observed_stat_examples_files/figure-html/unnamed-chunk-173-1.png new file mode 100644 index 00000000..f4d4797f Binary files /dev/null and b/docs/articles/observed_stat_examples_files/figure-html/unnamed-chunk-173-1.png differ diff --git a/docs/articles/observed_stat_examples_files/figure-html/unnamed-chunk-174-1.png b/docs/articles/observed_stat_examples_files/figure-html/unnamed-chunk-174-1.png new file mode 100644 index 00000000..0606b170 Binary files /dev/null and b/docs/articles/observed_stat_examples_files/figure-html/unnamed-chunk-174-1.png differ diff --git a/docs/articles/observed_stat_examples_files/figure-html/unnamed-chunk-179-1.png b/docs/articles/observed_stat_examples_files/figure-html/unnamed-chunk-179-1.png new file mode 100644 index 00000000..64a00b13 Binary files /dev/null and b/docs/articles/observed_stat_examples_files/figure-html/unnamed-chunk-179-1.png differ diff --git a/docs/articles/observed_stat_examples_files/figure-html/unnamed-chunk-19-1.png b/docs/articles/observed_stat_examples_files/figure-html/unnamed-chunk-19-1.png index 5909a038..9c6a85ca 100644 Binary files a/docs/articles/observed_stat_examples_files/figure-html/unnamed-chunk-19-1.png and b/docs/articles/observed_stat_examples_files/figure-html/unnamed-chunk-19-1.png differ diff --git a/docs/articles/observed_stat_examples_files/figure-html/unnamed-chunk-24-1.png b/docs/articles/observed_stat_examples_files/figure-html/unnamed-chunk-24-1.png new file mode 100644 index 00000000..7c186f97 Binary files /dev/null and b/docs/articles/observed_stat_examples_files/figure-html/unnamed-chunk-24-1.png differ diff --git a/docs/articles/observed_stat_examples_files/figure-html/unnamed-chunk-30-1.png b/docs/articles/observed_stat_examples_files/figure-html/unnamed-chunk-30-1.png index 5b4ffe85..e32337db 100644 Binary files a/docs/articles/observed_stat_examples_files/figure-html/unnamed-chunk-30-1.png and b/docs/articles/observed_stat_examples_files/figure-html/unnamed-chunk-30-1.png differ diff --git a/docs/articles/observed_stat_examples_files/figure-html/unnamed-chunk-35-1.png b/docs/articles/observed_stat_examples_files/figure-html/unnamed-chunk-35-1.png index cca8281a..7dc222c6 100644 Binary files a/docs/articles/observed_stat_examples_files/figure-html/unnamed-chunk-35-1.png and b/docs/articles/observed_stat_examples_files/figure-html/unnamed-chunk-35-1.png differ diff --git a/docs/articles/observed_stat_examples_files/figure-html/unnamed-chunk-40-1.png b/docs/articles/observed_stat_examples_files/figure-html/unnamed-chunk-40-1.png index 2e4b1f5f..dc28f3ce 100644 Binary files a/docs/articles/observed_stat_examples_files/figure-html/unnamed-chunk-40-1.png and b/docs/articles/observed_stat_examples_files/figure-html/unnamed-chunk-40-1.png differ diff --git a/docs/articles/observed_stat_examples_files/figure-html/unnamed-chunk-44-1.png b/docs/articles/observed_stat_examples_files/figure-html/unnamed-chunk-44-1.png new file mode 100644 index 00000000..34097149 Binary files /dev/null and b/docs/articles/observed_stat_examples_files/figure-html/unnamed-chunk-44-1.png differ diff --git a/docs/articles/observed_stat_examples_files/figure-html/unnamed-chunk-5-1.png b/docs/articles/observed_stat_examples_files/figure-html/unnamed-chunk-5-1.png new file mode 100644 index 00000000..bef89d52 Binary files /dev/null and b/docs/articles/observed_stat_examples_files/figure-html/unnamed-chunk-5-1.png differ diff --git a/docs/articles/observed_stat_examples_files/figure-html/unnamed-chunk-50-1.png b/docs/articles/observed_stat_examples_files/figure-html/unnamed-chunk-50-1.png new file mode 100644 index 00000000..57d3d397 Binary files /dev/null and b/docs/articles/observed_stat_examples_files/figure-html/unnamed-chunk-50-1.png differ diff --git a/docs/articles/observed_stat_examples_files/figure-html/unnamed-chunk-51-1.png b/docs/articles/observed_stat_examples_files/figure-html/unnamed-chunk-51-1.png index 109d7b15..b738c5b7 100644 Binary files a/docs/articles/observed_stat_examples_files/figure-html/unnamed-chunk-51-1.png and b/docs/articles/observed_stat_examples_files/figure-html/unnamed-chunk-51-1.png differ diff --git a/docs/articles/observed_stat_examples_files/figure-html/unnamed-chunk-52-1.png b/docs/articles/observed_stat_examples_files/figure-html/unnamed-chunk-52-1.png new file mode 100644 index 00000000..35d36ff0 Binary files /dev/null and b/docs/articles/observed_stat_examples_files/figure-html/unnamed-chunk-52-1.png differ diff --git a/docs/articles/observed_stat_examples_files/figure-html/unnamed-chunk-58-1.png b/docs/articles/observed_stat_examples_files/figure-html/unnamed-chunk-58-1.png index f4499cf6..88abfbd7 100644 Binary files a/docs/articles/observed_stat_examples_files/figure-html/unnamed-chunk-58-1.png and b/docs/articles/observed_stat_examples_files/figure-html/unnamed-chunk-58-1.png differ diff --git a/docs/articles/observed_stat_examples_files/figure-html/unnamed-chunk-59-1.png b/docs/articles/observed_stat_examples_files/figure-html/unnamed-chunk-59-1.png new file mode 100644 index 00000000..aa6a8cde Binary files /dev/null and b/docs/articles/observed_stat_examples_files/figure-html/unnamed-chunk-59-1.png differ diff --git a/docs/articles/observed_stat_examples_files/figure-html/unnamed-chunk-60-1.png b/docs/articles/observed_stat_examples_files/figure-html/unnamed-chunk-60-1.png new file mode 100644 index 00000000..aa6a8cde Binary files /dev/null and b/docs/articles/observed_stat_examples_files/figure-html/unnamed-chunk-60-1.png differ diff --git a/docs/articles/observed_stat_examples_files/figure-html/unnamed-chunk-67-1.png b/docs/articles/observed_stat_examples_files/figure-html/unnamed-chunk-67-1.png new file mode 100644 index 00000000..691ff5ae Binary files /dev/null and b/docs/articles/observed_stat_examples_files/figure-html/unnamed-chunk-67-1.png differ diff --git a/docs/articles/observed_stat_examples_files/figure-html/unnamed-chunk-68-1.png b/docs/articles/observed_stat_examples_files/figure-html/unnamed-chunk-68-1.png index b6217bde..38045fa6 100644 Binary files a/docs/articles/observed_stat_examples_files/figure-html/unnamed-chunk-68-1.png and b/docs/articles/observed_stat_examples_files/figure-html/unnamed-chunk-68-1.png differ diff --git a/docs/articles/observed_stat_examples_files/figure-html/unnamed-chunk-69-1.png b/docs/articles/observed_stat_examples_files/figure-html/unnamed-chunk-69-1.png index 427b4df8..738bcc82 100644 Binary files a/docs/articles/observed_stat_examples_files/figure-html/unnamed-chunk-69-1.png and b/docs/articles/observed_stat_examples_files/figure-html/unnamed-chunk-69-1.png differ diff --git a/docs/articles/observed_stat_examples_files/figure-html/unnamed-chunk-75-1.png b/docs/articles/observed_stat_examples_files/figure-html/unnamed-chunk-75-1.png new file mode 100644 index 00000000..f54dc748 Binary files /dev/null and b/docs/articles/observed_stat_examples_files/figure-html/unnamed-chunk-75-1.png differ diff --git a/docs/articles/observed_stat_examples_files/figure-html/unnamed-chunk-81-1.png b/docs/articles/observed_stat_examples_files/figure-html/unnamed-chunk-81-1.png index 5f4d0e54..46860cc0 100644 Binary files a/docs/articles/observed_stat_examples_files/figure-html/unnamed-chunk-81-1.png and b/docs/articles/observed_stat_examples_files/figure-html/unnamed-chunk-81-1.png differ diff --git a/docs/articles/observed_stat_examples_files/figure-html/unnamed-chunk-82-1.png b/docs/articles/observed_stat_examples_files/figure-html/unnamed-chunk-82-1.png new file mode 100644 index 00000000..f1aade7a Binary files /dev/null and b/docs/articles/observed_stat_examples_files/figure-html/unnamed-chunk-82-1.png differ diff --git a/docs/articles/observed_stat_examples_files/figure-html/unnamed-chunk-83-1.png b/docs/articles/observed_stat_examples_files/figure-html/unnamed-chunk-83-1.png new file mode 100644 index 00000000..daf09416 Binary files /dev/null and b/docs/articles/observed_stat_examples_files/figure-html/unnamed-chunk-83-1.png differ diff --git a/docs/articles/observed_stat_examples_files/figure-html/unnamed-chunk-88-1.png b/docs/articles/observed_stat_examples_files/figure-html/unnamed-chunk-88-1.png new file mode 100644 index 00000000..04952ec9 Binary files /dev/null and b/docs/articles/observed_stat_examples_files/figure-html/unnamed-chunk-88-1.png differ diff --git a/docs/articles/observed_stat_examples_files/figure-html/unnamed-chunk-94-1.png b/docs/articles/observed_stat_examples_files/figure-html/unnamed-chunk-94-1.png index abda9b8a..8469524b 100644 Binary files a/docs/articles/observed_stat_examples_files/figure-html/unnamed-chunk-94-1.png and b/docs/articles/observed_stat_examples_files/figure-html/unnamed-chunk-94-1.png differ diff --git a/docs/articles/observed_stat_examples_files/figure-html/unnamed-chunk-95-1.png b/docs/articles/observed_stat_examples_files/figure-html/unnamed-chunk-95-1.png index db71021b..e89c374c 100644 Binary files a/docs/articles/observed_stat_examples_files/figure-html/unnamed-chunk-95-1.png and b/docs/articles/observed_stat_examples_files/figure-html/unnamed-chunk-95-1.png differ diff --git a/docs/articles/observed_stat_examples_files/figure-html/unnamed-chunk-96-1.png b/docs/articles/observed_stat_examples_files/figure-html/unnamed-chunk-96-1.png new file mode 100644 index 00000000..36e26df1 Binary files /dev/null and b/docs/articles/observed_stat_examples_files/figure-html/unnamed-chunk-96-1.png differ diff --git a/docs/articles/observed_stat_examples_files/header-attrs-2.8/header-attrs.js b/docs/articles/observed_stat_examples_files/header-attrs-2.8/header-attrs.js new file mode 100644 index 00000000..dd57d92e --- /dev/null +++ b/docs/articles/observed_stat_examples_files/header-attrs-2.8/header-attrs.js @@ -0,0 +1,12 @@ +// Pandoc 2.9 adds attributes on both header and div. We remove the former (to +// be compatible with the behavior of Pandoc < 2.8). +document.addEventListener('DOMContentLoaded', function(e) { + var hs = document.querySelectorAll("div.section[class*='level'] > :first-child"); + var i, h, a; + for (i = 0; i < hs.length; i++) { + h = hs[i]; + if (!/^h[1-6]$/i.test(h.tagName)) continue; // it should be a header h1-h6 + a = h.attributes; + while (a.length > 0) h.removeAttribute(a[0].name); + } +}); diff --git a/docs/articles/t_test.html b/docs/articles/t_test.html index 4ce3c82a..42ca8409 100644 --- a/docs/articles/t_test.html +++ b/docs/articles/t_test.html @@ -12,21 +12,30 @@ - + - - - - + + + + + + + +
@@ -105,13 +109,13 @@ -
+
@@ -120,161 +124,184 @@

Tidy t-Tests with infer

-Introduction

-

In this vignette, we’ll walk through conducting t-tests using infer. We’ll start out with a 1-sample t-test, which compares a sample mean to a hypothesized true mean value. Then, we’ll discuss paired t-tests, which are a special use case of 1-sample t-tests, and evaluate whether differences in paired values (e.g. some measure taken of a person before and after an experiment) differ from 0. Finally, we’ll wrap up with 2-sample t-tests, testing the difference in means of two populations using a sample of data drawn from them.

+Introduction

+

In this vignette, we’ll walk through conducting \(t\)-tests and their randomization-based analogue using infer. We’ll start out with a 1-sample \(t\)-test, which compares a sample mean to a hypothesized true mean value. Then, we’ll discuss paired \(t\)-tests, which are a special use case of 1-sample \(t\)-tests, and evaluate whether differences in paired values (e.g. some measure taken of a person before and after an experiment) differ from 0. Finally, we’ll wrap up with 2-sample \(t\)-tests, testing the difference in means of two populations using a sample of data drawn from them.

Throughout this vignette, we’ll make use of the gss dataset supplied by infer, which contains a sample of data from the General Social Survey. See ?gss for more information on the variables included and their source. Note that this data (and our examples on it) are for demonstration purposes only, and will not necessarily provide accurate estimates unless weighted properly. For these examples, let’s suppose that this dataset is a representative sample of a population we want to learn about: American adults. The data looks like this:

-
dplyr::glimpse(gss)
+
+dplyr::glimpse(gss)
## Rows: 500
 ## Columns: 11
-## $ year    <dbl> 2014, 1994, 1998, 1996, 1994, 1996, 1990, 2016, 2000, 1998, 2…
-## $ age     <dbl> 36, 34, 24, 42, 31, 32, 48, 36, 30, 33, 21, 30, 38, 49, 25, 5…
-## $ sex     <fct> male, female, male, male, male, female, female, female, femal…
-## $ college <fct> degree, no degree, degree, no degree, degree, no degree, no d…
-## $ partyid <fct> ind, rep, ind, ind, rep, rep, dem, ind, rep, dem, dem, ind, d…
-## $ hompop  <dbl> 3, 4, 1, 4, 2, 4, 2, 1, 5, 2, 4, 3, 4, 4, 2, 2, 3, 2, 1, 2, 5…
-## $ hours   <dbl> 50, 31, 40, 40, 40, 53, 32, 20, 40, 40, 23, 52, 38, 72, 48, 4…
-## $ income  <ord> $25000 or more, $20000 - 24999, $25000 or more, $25000 or mor…
-## $ class   <fct> middle class, working class, working class, working class, mi…
-## $ finrela <fct> below average, below average, below average, above average, a…
-## $ weight  <dbl> 0.8960, 1.0825, 0.5501, 1.0864, 1.0825, 1.0864, 1.0627, 0.478…
+## $ year <dbl> 2014, 1994, 1998, 1996, 1994, 1996, 1990, 2016, 2000, 1998, 20… +## $ age <dbl> 36, 34, 24, 42, 31, 32, 48, 36, 30, 33, 21, 30, 38, 49, 25, 56… +## $ sex <fct> male, female, male, male, male, female, female, female, female… +## $ college <fct> degree, no degree, degree, no degree, degree, no degree, no de… +## $ partyid <fct> ind, rep, ind, ind, rep, rep, dem, ind, rep, dem, dem, ind, de… +## $ hompop <dbl> 3, 4, 1, 4, 2, 4, 2, 1, 5, 2, 4, 3, 4, 4, 2, 2, 3, 2, 1, 2, 5,… +## $ hours <dbl> 50, 31, 40, 40, 40, 53, 32, 20, 40, 40, 23, 52, 38, 72, 48, 40… +## $ income <ord> $25000 or more, $20000 - 24999, $25000 or more, $25000 or more… +## $ class <fct> middle class, working class, working class, working class, mid… +## $ finrela <fct> below average, below average, below average, above average, ab… +## $ weight <dbl> 0.8960, 1.0825, 0.5501, 1.0864, 1.0825, 1.0864, 1.0627, 0.4785…

-1-Sample t-Test

-

The 1-sample t-test can be used to test whether a sample of continuous data could have plausibly come from a population with a specified mean. As an example, we’ll test whether the average American adult works 40 hours a week using data from the gss. To do so, we make use of the hours variable, giving the number of hours that respondents reported having worked in the previous week. The distribution of hours in the observed data looks like this:

+1-Sample t-Test +

The 1-sample \(t\)-test can be used to test whether a sample of continuous data could have plausibly come from a population with a specified mean.

+

As an example, we’ll test whether the average American adult works 40 hours a week using data from the gss. To do so, we make use of the hours variable, giving the number of hours that respondents reported having worked in the previous week. The distribution of hours in the observed data looks like this:

-

Note the warning about missing values—many respondents’ entries are missing. If we were actually carrying out this hypothesis test, we might look further into how this data was collected; it’s possible that some of the missing values should actually be 0 hours.

-

In general, though, it looks like most respondents reported having worked 40 hours, but there’s quite a bit of variability. Let’s test whether we have evidence that the true mean number of hours that Americans work per week is 40.

+

It looks like most respondents reported having worked 40 hours, but there’s quite a bit of variability. Let’s test whether we have evidence that the true mean number of hours that Americans work per week is 40.

+

infer’s randomization-based analogue to the 1-sample \(t\)-test is a 1-sample mean test. We’ll start off showcasing that test before demonstrating how to carry out a theory-based \(t\)-test with the package.

First, to calculate the observed statistic, we can use specify() and calculate().

-
# calculate the observed statistic
-observed_statistic <- gss %>%
-  specify(response = hours) %>%
-  hypothesize(null = "point", mu = 40) %>%
-  calculate(stat = "t")
-

The observed statistic is 2.0852. Now, we want to compare this statistic to a null distribution, generated under the assumption that the mean was actually 40, to get a sense of how likely it would be for us to see this observed statistic if the true number of hours worked per week in the population was really 40.

+
+# calculate the observed statistic
+observed_statistic <- gss %>%
+  specify(response = hours) %>%
+  calculate(stat = "mean")
+

The observed statistic is 41.382. Now, we want to compare this statistic to a null distribution, generated under the assumption that the mean was actually 40, to get a sense of how likely it would be for us to see this observed mean if the true number of hours worked per week in the population was really 40.

We can generate the null distribution using the bootstrap. In the bootstrap, for each replicate, a sample of size equal to the input sample size is drawn (with replacement) from the input sample data. This allows us to get a sense of how much variability we’d expect to see in the entire population so that we can then understand how unlikely our sample mean would be.

-
# generate the null distribution
-null_distribution_1_sample <- gss %>%
-  specify(response = hours) %>%
-  hypothesize(null = "point", mu = 40) %>%
-  generate(reps = 1000, type = "bootstrap") %>%
-  calculate(stat = "t")
+
+# generate the null distribution
+null_dist_1_sample <- gss %>%
+  specify(response = hours) %>%
+  hypothesize(null = "point", mu = 40) %>%
+  generate(reps = 1000, type = "bootstrap") %>%
+  calculate(stat = "mean")

To get a sense for what these distributions look like, and where our observed statistic falls, we can use visualize():

-
# visualize the null distribution and test statistic!
-null_distribution_1_sample %>%
-  visualize() +
-  shade_p_value(observed_statistic,
-                direction = "two-sided")
+
+# visualize the null distribution and test statistic!
+null_dist_1_sample %>%
+  visualize() + 
+  shade_p_value(observed_statistic,
+                direction = "two-sided")

-

It looks like our observed statistic of 2.0852 would be relatively unlikely if the true mean was actually 40 hours a week. More exactly, we can calculate the p-value:

-
# calculate the p value from the test statistic and null distribution
-p_value_1_sample <- null_distribution_1_sample %>%
-  get_p_value(obs_stat = observed_statistic,
-              direction = "two-sided")
+

It looks like our observed mean of 41.382 would be relatively unlikely if the true mean was actually 40 hours a week. More exactly, we can calculate the p-value:

+
+# calculate the p value from the test statistic and null distribution
+p_value_1_sample <- null_dist_1_sample %>%
+  get_p_value(obs_stat = observed_statistic,
+              direction = "two-sided")
 
-p_value_1_sample
-
## # A tibble: 1 x 1
+p_value_1_sample
+
## # A tibble: 1 × 1
 ##   p_value
 ##     <dbl>
-## 1   0.032
-

Thus, if the true mean number of hours worked per week was really 40, the probability that we would see a test statistic as or more extreme than 2.0852 is approximately 0.032.

-

Note that, equivalently to the steps shown above, the package supplies a wrapper function, t_test, to carry out 1-sample t-Tests on tidy data. The syntax looks like this:

-
t_test(gss, response = hours, mu = 40)
-
## # A tibble: 1 x 6
-##   statistic  t_df p_value alternative lower_ci upper_ci
-##       <dbl> <dbl>   <dbl> <chr>          <dbl>    <dbl>
-## 1      2.09   499  0.0376 two.sided       40.1     42.7
-

Also, to just calculate the observed statistic, the package supplies the wrapper t_stat().

-
t_stat(gss, response = hours, mu = 40)
-
##     t 
-## 2.085
+## 1 0.026
+

Thus, if the true mean number of hours worked per week was really 40, our approximation of the probability that we would see a test statistic as or more extreme than 41.382 is approximately 0.026.

+

Analogously to the steps shown above, the package supplies a wrapper function, t_test, to carry out 1-sample \(t\)-tests on tidy data. Rather than using randomization, the wrappers carry out the theory-based \(t\)-test. The syntax looks like this:

+
+t_test(gss, response = hours, mu = 40)
+
## # A tibble: 1 × 7
+##   statistic  t_df p_value alternative estimate lower_ci upper_ci
+##       <dbl> <dbl>   <dbl> <chr>          <dbl>    <dbl>    <dbl>
+## 1      2.09   499  0.0376 two.sided       41.4     40.1     42.7
+

An alternative approach to the t_test() wrapper is to calculate the observed statistic with an infer pipeline and then supply it to the pt function from base R.

+
+# calculate the observed statistic
+observed_statistic <- gss %>%
+  specify(response = hours) %>%
+  hypothesize(null = "point", mu = 40) %>%
+  calculate(stat = "t") %>%
+  dplyr::pull()
+

Note that this pipeline to calculate an observed statistic includes a call to hypothesize() since the \(t\) statistic requires a hypothesized mean value.

+

Then, juxtaposing that \(t\) statistic with its associated distribution using the pt function:

+
+pt(observed_statistic, df = nrow(gss) - 1, lower.tail = FALSE)*2
+
##       t 
+## 0.03756
+

Note that the resulting \(t\)-statistics from these two theory-based approaches are the same.

-Paired t-Test

-

You might be interested in running a paired t-Test. Paired t-tests can be used in situations when there is a natural pairing between values in distributions—a common example would be two columns, before and after, say, that contain measurements from a patient before and after some treatment. To compare these two distributions, then, we’re not necessarily interested in how the two distributions look different altogether, but how these two measurements from each individual change across time. (Pairings don’t necessarily have to be over time; another common usage is measurements from two married people, for example.) Thus, we can create a new column (see mutate() from the dplyr package if you’re not sure how to do this) that is the difference between the two: difference = after - before, and then examine this distribution to see how each individuals’ measurements changed over time.

-

Once we’ve mutate()d that new difference column, we can run a 1-sample t-Test on it, where our null hypothesis is that mu = 0 (i.e. the difference between these measurements before and after treatment is, on average, 0). To do so, we’d use the procedure outlined in the above section.

+Paired t-Test +

You might be interested in running a paired \(t\)-test. Paired \(t\)-tests can be used in situations when there is a natural pairing between values in distributions—a common example would be two columns, before and after, say, that contain measurements from a patient before and after some treatment. To compare these two distributions, then, we’re not necessarily interested in how the two distributions look different altogether, but how these two measurements from each individual change across time. (Pairings don’t necessarily have to be over time; another common usage is measurements from two married people, for example.) Thus, we can create a new column (see mutate() from the dplyr package if you’re not sure how to do this) that is the difference between the two: difference = after - before, and then examine this distribution to see how each individuals’ measurements changed over time.

+

Once we’ve mutate()d that new difference column, we can run a 1-sample \(t\)-test on it, where our null hypothesis is that mu = 0 (i.e. the difference between these measurements before and after treatment is, on average, 0). To do so, we’d use the procedure outlined in the above section.

-2-Sample t-Test

-

2-Sample t-Tests evaluate the difference in mean values of two populations using data randomly-sampled from the population that approximately follows a normal distribution. As an example, we’ll test if Americans work the same number of hours a week regardless of whether they have a college degree or not using data from the gss. The college and hours variables allow us to do so:

+2-Sample t-Test +

2-Sample \(t\)-tests evaluate the difference in mean values of two populations using data randomly-sampled from the population that approximately follows a normal distribution. As an example, we’ll test if Americans work the same number of hours a week regardless of whether they have a college degree or not using data from the gss. The college and hours variables allow us to do so:

It looks like both of these distributions are centered near 40 hours a week, but the distribution for those with a degree is slightly right skewed.

Again, note the warning about missing values—many respondents’ values are missing. If we were actually carrying out this hypothesis test, we might look further into how this data was collected; it’s possible that whether or not a value in either of these columns is missing is related to what that value would be.

-

Now, to calculate the observed \(t\) statistic, we can use specify() and calculate().

-
# calculate the observed statistic
-observed_statistic <- gss %>%
-  specify(hours ~ college) %>%
-  calculate(stat = "t", order = c("degree", "no degree"))
+

infer’s randomization-based analogue to the 2-sample \(t\)-test is a difference in means test. We’ll start off showcasing that test before demonstrating how to carry out a theory-based \(t\)-test with the package.

+

As with the one-sample test, to calculate the observed difference in means, we can use specify() and calculate().

+
+# calculate the observed statistic
+observed_statistic <- gss %>%
+  specify(hours ~ college) %>%
+  calculate(stat = "diff in means", order = c("degree", "no degree"))
 
-observed_statistic
-
## # A tibble: 1 x 1
+observed_statistic
+
## Response: hours (numeric)
+## Explanatory: college (factor)
+## # A tibble: 1 × 1
 ##    stat
 ##   <dbl>
-## 1  1.12
-

Note that, in the line specify(hours ~ college), we could have swapped this out with the syntax specify(response = hours, explanatory = college)!

+## 1 1.54
+

Note that, in the line specify(hours ~ college), we could have swapped this out with the syntax specify(response = hours, explanatory = college)!

The order argument in that calculate line gives the order to subtract the mean values in: in our case, we’re taking the mean number of hours worked by those with a degree minus the mean number of hours worked by those without a degree; a positive difference, then, would mean that people with degrees worked more than those without a degree.

-

Now, we want to compare this \(t\)-statistic to a null distribution, generated under the assumption that the number of hours worked a week has no relationship with whether or not one has a college degree, to get a sense of how likely it would be for us to see this \(t\)-statistic if there were really no relationship between these two variables.

+

Now, we want to compare this difference in means to a null distribution, generated under the assumption that the number of hours worked a week has no relationship with whether or not one has a college degree, to get a sense of how likely it would be for us to see this observed difference in means if there were really no relationship between these two variables.

We can generate the null distribution using permutation, where, for each replicate, each value of degree status will be randomly reassigned (without replacement) to a new number of hours worked per week in the sample in order to break any association between the two.

-
# generate the null distribution with randomization
-null_distribution_2_sample_permute <- gss %>%
-  specify(hours ~ college) %>%
-  hypothesize(null = "independence") %>%
-  generate(reps = 1000, type = "permute") %>%
-  calculate(stat = "t", order = c("degree", "no degree"))
-

We can also generate the null distribution using the theoretical \(t\) distribution. The code to do that looks like this:

-
# generate the null distribution with the theoretical t
-null_distribution_2_sample_theoretical <- gss %>%
-  specify(hours ~ college) %>%
-  hypothesize(null = "independence") %>%
-  # generate() isn't used for the theoretical version!
-  calculate(stat = "t", order = c("degree", "no degree"))
-

Again, note that, in the lines specify(hours ~ college) in the above two chunks, we could have used the syntax specify(response = hours, explanatory = college) instead!

-

To get a sense for what these distributions look like, and where our observed statistic falls, we can use visualize(). First, to do this with the randomization-based distribution:

-
# visualize the randomization-based null distribution and test statistic!
-null_distribution_2_sample_permute %>%
-  visualize() +
-  shade_p_value(observed_statistic,
-                direction = "two-sided")
+
+# generate the null distribution with randomization
+null_dist_2_sample <- gss %>%
+  specify(hours ~ college) %>%
+  hypothesize(null = "independence") %>%
+  generate(reps = 1000, type = "permute") %>%
+  calculate(stat = "diff in means", order = c("degree", "no degree"))
+

Again, note that, in the lines specify(hours ~ college) in the above chunk, we could have used the syntax specify(response = hours, explanatory = college) instead!

+

To get a sense for what these distributions look like, and where our observed statistic falls, we can use visualize().

+
+# visualize the randomization-based null distribution and test statistic!
+null_dist_2_sample %>%
+  visualize() + 
+  shade_p_value(observed_statistic,
+                direction = "two-sided")

-

The syntax to do the same thing with the theory-based null distribution looks similar—we just pipe in the theory-based null distribution instead, and drop in method = "theoretical" to visualize().

-
# visualize the theoretical null distribution and test statistic!
-null_distribution_2_sample_theoretical %>%
-  visualize(method = "theoretical") +
-  shade_p_value(observed_statistic,
-                direction = "two-sided")
-

-

We could also visualize both null distributions to get a sense of how the randomization-based and theoretical versions relate. To do so, start off with the randomization-based null distribution, and then drop in method = "both" to visualize().

-
# visualize both null distributions and test statistic!
-null_distribution_2_sample_permute %>%
-  visualize(method = "both") +
-  shade_p_value(observed_statistic,
-                direction = "two-sided")
-

-

Regardless, it looks like our observed statistic of 1.1193 would be really unlikely if there was truly no relationship between degree status and number of hours worked. More exactly, we can calculate the p-value; theoretical p-values are not yet supported, so we’ll use the randomization-based null distribution to do calculate the p-value.

-
# calculate the p value from the randomization-based null 
+

It looks like our observed statistic of 1.5384 would be unlikely if there was truly no relationship between degree status and number of hours worked. More exactly, we can calculate the p-value; theoretical p-values are not yet supported, so we’ll use the randomization-based null distribution to do calculate the p-value.

+
+# calculate the p value from the randomization-based null 
 # distribution and the observed statistic
-p_value_2_sample <- null_distribution_2_sample_permute %>%
-  get_p_value(obs_stat = observed_statistic,
-              direction = "two-sided")
+p_value_2_sample <- null_dist_2_sample %>%
+  get_p_value(obs_stat = observed_statistic,
+              direction = "two-sided")
 
-p_value_2_sample
-
## # A tibble: 1 x 1
+p_value_2_sample
+
## # A tibble: 1 × 1
 ##   p_value
 ##     <dbl>
-## 1   0.268
-

Thus, if there were really no relationship between the number of hours worked a week and whether one has a college degree, the probability that we would see a statistic as or more extreme than 1.1193 is approximately 0.268.

-

Note that, equivalently to the steps shown above, the package supplies a wrapper function, t_test, to carry out 2-sample t-Tests on tidy data. The syntax looks like this:

-
t_test(x = gss,
-       formula = hours ~ college,
-       order = c("degree", "no degree"),
-       alternative = "two-sided")
-
## # A tibble: 1 x 6
-##   statistic  t_df p_value alternative lower_ci upper_ci
-##       <dbl> <dbl>   <dbl> <chr>          <dbl>    <dbl>
-## 1      1.12  366.   0.264 two.sided      -1.16     4.24
+## 1 0.3
+

Thus, if there were really no relationship between the number of hours worked a week and whether one has a college degree, the probability that we would see a statistic as or more extreme than 1.5384 is approximately 0.3.

+

Note that, similarly to the steps shown above, the package supplies a wrapper function, t_test, to carry out 2-sample \(t\)-tests on tidy data. The syntax looks like this:

+
+t_test(x = gss, 
+       formula = hours ~ college, 
+       order = c("degree", "no degree"),
+       alternative = "two-sided")
+
## # A tibble: 1 × 7
+##   statistic  t_df p_value alternative estimate lower_ci upper_ci
+##       <dbl> <dbl>   <dbl> <chr>          <dbl>    <dbl>    <dbl>
+## 1      1.12  366.   0.264 two.sided       1.54    -1.16     4.24

In the above example, we specified the relationship with the syntax formula = hours ~ college; we could have also written response = hours, explanatory = college.

+

An alternative approach to the t_test() wrapper is to calculate the observed statistic with an infer pipeline and then supply it to the pt function from base R. We can calculate the statistic as before, switching out the stat = "diff in means" argument with stat = "t".

+
+# calculate the observed statistic
+observed_statistic <- gss %>%
+  specify(hours ~ college) %>%
+  hypothesize(null = "point", mu = 40) %>%
+  calculate(stat = "t", order = c("degree", "no degree")) %>%
+  dplyr::pull()
+
+observed_statistic
+
##     t 
+## 1.119
+

Note that this pipeline to calculate an observed statistic includes hypothesize() since the \(t\) statistic requires a hypothesized mean value.

+

Then, juxtaposing that \(t\) statistic with its associated distribution using the pt function:

+
+pt(observed_statistic, df = nrow(gss) - 2, lower.tail = FALSE)*2
+
##      t 
+## 0.2635
+

Note that the results from these two theory-based approaches are nearly the same.

@@ -286,32 +313,24 @@

-

@@ -155,34 +161,34 @@ - + + @@ -222,36 +229,24 @@

Authors

-
- - + diff --git a/docs/index.html b/docs/index.html index 2a53b0f9..9d95391a 100644 --- a/docs/index.html +++ b/docs/index.html @@ -12,22 +12,31 @@ - + - - - - + + + + + + + +
@@ -110,10 +114,13 @@
-
- + + + + +

The objective of this package is to perform statistical inference using an expressive statistical grammar that coheres with the tidyverse design framework. The package is centered around 4 main verbs, supplemented with many utilities to visualize and extract value from their outputs.

  • @@ -126,72 +133,83 @@ calculate() allows you to calculate a distribution of statistics from the generated data to form the null distribution.

To learn more about the principles underlying the package design, see vignette("infer").

-

+
+A diagram showing four steps to carry out randomization-based inference: specify hypothesis, generate data, calculate statistic, and visualize. From left to right, each step is connected by an arrow, while the diagram indicates that generating data and calculating statistics can happen iteratively.

+

+
+

If you’re interested in learning more about randomization-based statistical inference generally, including applied examples of this package, we recommend checking out Statistical Inference Via Data Science: A ModernDive Into R and the Tidyverse and Introduction to Modern Statistics.

-Installation

+Installation

To install the current stable version of infer from CRAN:

- +

To install the developmental stable version of infer, make sure to install remotes first. The pkgdown website for this version is at infer.tidymodels.org.

-
install.packages("remotes")
-remotes::install_github("tidymodels/infer")
-

To install the experimental version of infer (do so at your own risk), make sure to install remotes first.

-
install.packages("remotes")
-remotes::install_github("tidymodels/infer", ref = "develop")
-

To see documentation for the experimental version of infer, the pkgdown site is available at https://infer-dev.netlify.com.

+
+install.packages("remotes")
+remotes::install_github("tidymodels/infer")

-Contributing

+Contributing
-

We welcome others helping us make this package as user-friendly and efficient as possible. Please review our contributing and conduct guidelines. Please see the open issues for more specific fixes/additions needing attention. By participating in this project you agree to abide by its terms.

+

We welcome others helping us make this package as user-friendly and efficient as possible. Please review our contributing and conduct guidelines. By participating in this project you agree to abide by its terms.

+

For questions and discussions about tidymodels packages, modeling, and machine learning, please post on RStudio Community. If you think you have encountered a bug, please submit an issue. Either way, learn how to create and share a reprex (a minimal, reproducible example), to clearly communicate about your code. Check out further details on contributing guidelines for tidymodels packages and how to get help.

-Examples

+Examples
-

These examples are pulled from the “Full infer Pipeline Examples” vignette, accessible by calling vignette("observed_stat_examples"). They make use of the gss dataset supplied by the package, providing a sample of data from the General Social Survey. The data looks like this:

-
# load in the dataset
-data(gss)
+

These examples are pulled from the “Full infer Pipeline Examples” vignette, accessible by calling vignette("observed_stat_examples"). They make use of the gss dataset supplied by the package, providing a sample of data from the General Social Survey. The data looks like this:

+
+# load in the dataset
+data(gss)
 
 # take a glimpse at it
-str(gss)
-
## Classes 'tbl_df', 'tbl' and 'data.frame':    500 obs. of  11 variables:
-##  $ year   : num  2014 1994 1998 1996 1994 ...
-##  $ age    : num  36 34 24 42 31 32 48 36 30 33 ...
-##  $ sex    : Factor w/ 2 levels "male","female": 1 2 1 1 1 2 2 2 2 2 ...
-##  $ college: Factor w/ 2 levels "no degree","degree": 2 1 2 1 2 1 1 2 2 1 ...
-##  $ partyid: Factor w/ 5 levels "dem","ind","rep",..: 2 3 2 2 3 3 1 2 3 1 ...
-##  $ hompop : num  3 4 1 4 2 4 2 1 5 2 ...
-##  $ hours  : num  50 31 40 40 40 53 32 20 40 40 ...
-##  $ income : Ord.factor w/ 12 levels "lt $1000"<"$1000 to 2999"<..: 12 11 12 12 12 12 12 12 12 10 ...
-##  $ class  : Factor w/ 6 levels "lower class",..: 3 2 2 2 3 3 2 3 3 2 ...
-##  $ finrela: Factor w/ 6 levels "far below average",..: 2 2 2 4 4 3 2 4 3 1 ...
-##  $ weight : num  0.896 1.083 0.55 1.086 1.083 ...
+str(gss)
+
## tibble [500 × 11] (S3: tbl_df/tbl/data.frame)
+##  $ year   : num [1:500] 2014 1994 1998 1996 1994 ...
+##  $ age    : num [1:500] 36 34 24 42 31 32 48 36 30 33 ...
+##  $ sex    : Factor w/ 2 levels "male","female": 1 2 1 1 1 2 2 2 2 2 ...
+##  $ college: Factor w/ 2 levels "no degree","degree": 2 1 2 1 2 1 1 2 2 1 ...
+##  $ partyid: Factor w/ 5 levels "dem","ind","rep",..: 2 3 2 2 3 3 1 2 3 1 ...
+##  $ hompop : num [1:500] 3 4 1 4 2 4 2 1 5 2 ...
+##  $ hours  : num [1:500] 50 31 40 40 40 53 32 20 40 40 ...
+##  $ income : Ord.factor w/ 12 levels "lt $1000"<"$1000 to 2999"<..: 12 11 12 12 12 12 12 12 12 10 ...
+##  $ class  : Factor w/ 6 levels "lower class",..: 3 2 2 2 3 3 2 3 3 2 ...
+##  $ finrela: Factor w/ 6 levels "far below average",..: 2 2 2 4 4 3 2 4 3 1 ...
+##  $ weight : num [1:500] 0.896 1.083 0.55 1.086 1.083 ...

As an example, we’ll run an analysis of variance on age and partyid, testing whether the age of a respondent is independent of their political party affiliation.

Calculating the observed statistic,

-
F_hat <- gss %>%
-  specify(age ~ partyid) %>%
-  calculate(stat = "F")
+
+F_hat <- gss %>% 
+  specify(age ~ partyid) %>%
+  calculate(stat = "F")

Then, generating the null distribution,

-
null_distn <- gss %>%
-   specify(age ~ partyid) %>%
-   hypothesize(null = "independence") %>%
-   generate(reps = 1000, type = "permute") %>%
-   calculate(stat = "F")
+
+null_dist <- gss %>%
+   specify(age ~ partyid) %>%
+   hypothesize(null = "independence") %>%
+   generate(reps = 1000, type = "permute") %>%
+   calculate(stat = "F")

Visualizing the observed statistic alongside the null distribution,

-
visualize(null_distn) +
-  shade_p_value(obs_stat = F_hat, direction = "greater")
-

+
+visualize(null_dist) +
+  shade_p_value(obs_stat = F_hat, direction = "greater")
+
+A histogram showing a distribution of F statistics, right-tailed and centered around one. The x axis ranges from zero to five. The region of the histogram to the right of the observed statistic, just above two, is shaded red to represent the p-value.

+

+

Calculating the p-value from the null distribution and observed statistic,

-
null_distn %>%
-  get_p_value(obs_stat = F_hat, direction = "greater")
-
## # A tibble: 1 x 1
-##   p_value
-##     <dbl>
-## 1   0.056
-

Note that the formula and non-formula interfaces (i.e. age ~ partyid vs. response = age, explanatory = partyid) work for all implemented inference procedures in infer. Use whatever is more natural for you. If you will be doing modeling using functions like lm() and glm(), though, we recommend you begin to use the formula y ~ x notation as soon as possible.

+
+null_dist %>%
+  get_p_value(obs_stat = F_hat, direction = "greater")
+
## # A tibble: 1 x 1
+##   p_value
+##     <dbl>
+## 1   0.055
+

Note that the formula and non-formula interfaces (i.e. age ~ partyid vs. response = age, explanatory = partyid) work for all implemented inference procedures in infer. Use whatever is more natural for you. If you will be doing modeling using functions like lm() and glm(), though, we recommend you begin to use the formula y ~ x notation as soon as possible.

Other resources are available in the package vignettes! See vignette("observed_stat_examples") for more examples like the one above, and vignette("infer") for discussion of the underlying principles of the package design.

@@ -199,7 +217,7 @@

-