Practice Set 7.5

library(tidyverse)

Warning: package 'tidyverse' was built under R version 4.3.3

Warning: package 'tibble' was built under R version 4.3.3

Warning: package 'tidyr' was built under R version 4.3.3

Warning: package 'readr' was built under R version 4.3.3

Warning: package 'purrr' was built under R version 4.3.3

Warning: package 'dplyr' was built under R version 4.3.3

Warning: package 'stringr' was built under R version 4.3.3

Warning: package 'forcats' was built under R version 4.3.3

Warning: package 'lubridate' was built under R version 4.3.3

── Attaching core tidyverse packages ──────────────────────── tidyverse 2.0.0 ──
✔ dplyr     1.1.4     ✔ readr     2.1.5
✔ forcats   1.0.0     ✔ stringr   1.5.1
✔ ggplot2   4.0.1     ✔ tibble    3.2.1
✔ lubridate 1.9.4     ✔ tidyr     1.3.1
✔ purrr     1.0.2     
── Conflicts ────────────────────────────────────────── tidyverse_conflicts() ──
✖ dplyr::filter() masks stats::filter()
✖ dplyr::lag()    masks stats::lag()
ℹ Use the conflicted package (<http://conflicted.r-lib.org/>) to force all conflicts to become errors

library(purrrfect)

Attaching package: 'purrrfect'

The following objects are masked from 'package:base':

    replicate, tabulate

library(ggh4x)

Practice Set 7.5

Problem 4

Setting up the dataframe

q4df <- (parameters(~n,~mu,~sigma,c(4,8,15),c(-2,0,2),c(1,2,3))
    %>%add_trials(10000)
    %>%mutate(Y_n=pmap(list(n,mu,sigma),\(x,m,y) rnorm(x,m,y^2)),Y_hat = map_dbl(Y_n,\(y) mean(y)),S_n=pmap_dbl(list(n,Y_n,Y_hat),\(x,y,z) (sum((y-z)^2)/(x-1))^(1/2)))
    %>%mutate(T=pmap_dbl(list(n,mu,Y_hat,S_n),\(n,m,y,s) (y-m)/(s/(n)^(1/2))),ft = pmap_dbl(list(T,n),\(x,n) dt(x,n-1)))
)
q4df

# A tibble: 270,000 × 9
       n    mu sigma .trial Y_n       Y_hat   S_n       T     ft
   <dbl> <dbl> <dbl>  <dbl> <list>    <dbl> <dbl>   <dbl>  <dbl>
 1     4    -2     1      1 <dbl [4]> -1.58 1.06   0.791  0.252 
 2     4    -2     1      2 <dbl [4]> -2.07 0.630 -0.226  0.355 
 3     4    -2     1      3 <dbl [4]> -1.65 0.369  1.89   0.0765
 4     4    -2     1      4 <dbl [4]> -1.48 1.34   0.771  0.256 
 5     4    -2     1      5 <dbl [4]> -2.25 0.795 -0.636  0.285 
 6     4    -2     1      6 <dbl [4]> -2.11 1.23  -0.180  0.360 
 7     4    -2     1      7 <dbl [4]> -2.34 1.15  -0.592  0.295 
 8     4    -2     1      8 <dbl [4]> -1.92 0.622  0.271  0.350 
 9     4    -2     1      9 <dbl [4]> -2.33 0.499 -1.33   0.146 
10     4    -2     1     10 <dbl [4]> -1.98 0.880  0.0408 0.367 
# ℹ 269,990 more rows

Plotting analytic versus simulated variance

(ggplot(data = q4df) 
+ geom_histogram(aes(x = T,y = after_stat(density)),fill = "red",binwidth = .2,center = .01)
+ geom_line(aes(x = T,y = ft),color = "blue")
+ facet_nested(n~sigma+mu, labeller = label_both,scale = "free_y")
+ xlim(c(-5,5))
)

Warning: Removed 1479 rows containing non-finite outside the scale range
(`stat_bin()`).

Warning: Removed 54 rows containing missing values or values outside the scale range
(`geom_bar()`).

Warning: Removed 1479 rows containing missing values or values outside the scale range
(`geom_line()`).