Practice Set 7.4

library(purrrfect)

Attaching package: 'purrrfect'

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

    replicate, tabulate

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

Practice Set 7.4

Problem 5

Creating the dataset

q5df <- (parameters(~n,~sig,c(4,8,15),c(1,2,3))
  %>%add_trials(10000)
  %>%mutate(Y_n=pmap(list(n,sig),\(x,y) rnorm(x,0,y)),Y_hat = map(Y_n,\(y) mean(y)),S_n=pmap_dbl(list(n,Y_n,Y_hat),\(x,y,z) sum((y-z)^2)/(x-1)),S2= pmap_dbl(list(S_n,n,sig),\(x,y,z) ((y-1)*x)/(z^2)),Sf = pmap_dbl(list(S2,n),\(x,y) dgamma(x,(y-1)/2,1/2)))
)

Plotting Simulated versus Analytic PDF

(ggplot(data = q5df)
+ geom_histogram(aes(x = S2,y = after_stat(density)), fill = "gold",binwidth = .2,center = .01)
+ geom_line(aes(x=S2,y=Sf),color = "blue")
+ facet_grid(~n~sig,labeller =label_both,scales = "free")
+ theme_classic()
)

Table of simulated versus analytic expectations

(q5df %>% mutate(AE = pmap_dbl(list(n,sig),\(x,y) y^2))
    %>%summarize("E(S2)" = mean(S_n),"AE(S2)" = mean(AE),.by = c(n,sig))
)

# A tibble: 9 × 4
      n   sig `E(S2)` `AE(S2)`
  <dbl> <dbl>   <dbl>    <dbl>
1     4     1   1.02         1
2     4     2   3.99         4
3     4     3   8.95         9
4     8     1   1.00         1
5     8     2   4.05         4
6     8     3   8.99         9
7    15     1   0.995        1
8    15     2   3.99         4
9    15     3   9.00         9

table of simulated versus analytic variance’s

(q5df %>% mutate(AV = pmap_dbl(list(n,sig),\(x,y) 2*(y^4)/(x-1)))
    %>%summarize("E(S2)" = var(S_n),"AV(S2)" = mean(AV),.by = c(n,sig))
)

# A tibble: 9 × 4
      n   sig `E(S2)` `AV(S2)`
  <dbl> <dbl>   <dbl>    <dbl>
1     4     1   0.700    0.667
2     4     2  10.8     10.7  
3     4     3  53.6     54    
4     8     1   0.291    0.286
5     8     2   4.67     4.57 
6     8     3  23.1     23.1  
7    15     1   0.144    0.143
8    15     2   2.28     2.29 
9    15     3  11.5     11.6  

Problem 6

Creating dataset from exponential sample

q6df <- (parameters(~n,~lambda,c(10,20,30),c(.5,1,1.5))
  %>%add_trials(1000)
  %>%mutate(Y_n=pmap(list(n,lambda),\(x,y) rexp(x,y)),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)))
)

Plotting all combinations of\(S^2\) versus \(\overline{Y}\)

(ggplot(aes(x = S_n,y = Y_hat),data = q6df)
+ geom_point()
+ facet_grid(~n~lambda,labeller =label_both,scales = "free")
+ theme_classic()
)

As shown in above plots, \(S^2\) and \(\overline{Y}\) appear to not be independent as the points appear to be centered around a curved pattern for all combinations.