lab:projects:02robust_growth_curve_modeling_using_student-t_distribution:r_codes_related_to_real_data_analysis
## coda analysis
student<-read.coda('coda-uniform.txt','index-uniform.txt')
geweke.diag(student,.2,.5)
HPDinterval(student)
## history plot
## plot coda
coda<-read.table('/Volumes/~zzhang4/Private/zzy/research/Bayesian analysis of longitudinal data/Student t distribution/real-data-analysis/latent basis coda chain.txt')
par<-matrix(c(coda[,2]),ncol=11)
## history plot
par(mfrow=c(5,2), mar=c(3, 5, 1, 1) + 0.1)
plot.ts(par[,1], ylab='A[2]')
plot.ts(par[,2], ylab='A[3]')
plot.ts(par[,3], ylab='A[4]')
plot.ts(par[,9], ylab=quote(beta[L]))
plot.ts(par[,10], ylab=quote(beta[S]))
plot.ts(par[,4], ylab=quote(sigma[L]^2))
plot.ts(par[,5], ylab=quote(sigma[S]^2))
plot.ts(par[,6], ylab=quote(sigma[LS]))
plot.ts(par[,8], ylab=quote(phi))
plot.ts(par[,11], ylab=quote(k))
normal<-read.coda('coda-norm.txt','index-norm.txt')
geweke.diag(normal, .2,.5)
HPDinterval(normal)
id<-rep(1:7, each=38001)
x<-rep(2000:40000, 7)
par(mfrow=c(7,1), mar=c(3, 5, 1, 1) + 0.1)
plot(x[id==6], normal[id==6], type='l', ylab=quote(beta[L]), xlab='')
plot(x[id==7], normal[id==7], type='l', ylab=quote(beta[S]), xlab='')
plot(x[id==1], normal[id==1], type='l', ylab=quote(sigma[L]^2), xlab='')
plot(x[id==4], normal[id==4], type='l', ylab=quote(sigma[S]^2), xlab='')
plot(x[id==2], normal[id==2], type='l', ylab=quote(sigma[LS]), xlab='')
plot(x[id==5], normal[id==5], type='l', ylab=quote(phi), xlab='')
#####################
## Boxplots
y = matrix( c(
NA, 65.0, 50.0, 61.0, NA,
48.0, 44.0, 50.0, NA, 62.0,
60.0, 37.0, 42.0, 52.0, 53.0,
61.0, 85.0, 69.0, 71.0, 71.0,
NA, 38.0, 72.0, 69.0, 76.0,
59.0, 79.0, 61.0, 57.0, 70.0,
52.0, 72.0, 79.0, 69.0, NA,
73.0, 79.0, 82.0, 73.0, 90.0,
41.0, 67.0, NA, NA, 70.0,
69.0, 78.0, 72.0, 70.0, 74.0,
51.0, 45.0, 50.0, 70.0, 69.0,
63.0, 72.0, 69.0, 74.0, NA,
80.0, 85.0, 89.0, 92.0, 97.0,
78.0, 81.0, 96.0, 97.0, 81.0,
74.0, 50.0, NA, 62.0, 38.0,
80.0, 65.0, 85.0, 86.0, 86.0,
63.0, 72.0, 64.0, NA, 69.0,
68.0, 86.0, 73.0, 90.0, 95.0,
82.0, 92.0, 95.0, 97.0, 98.0,
61.0, 58.0, 79.0, 72.0, 72.0,
76.0, 55.0, 69.0, 82.0, 58.0,
80.0, 81.0, 87.0, 85.0, 92.0,
50.0, 48.0, 61.0, 73.0, 50.0,
64.0, NA, 57.0, 65.0, 55.0,
93.0, 93.0, 95.0, 97.0, 100.0,
50.0, 58.0, 77.0, 76.0, 58.0,
49.0, 54.0, 72.0, 65.0, 64.0,
58.0, 62.0, 81.0, 81.0, 84.0,
60.0, 70.0, 72.0, 82.0, 83.0,
67.0, 90.0, 84.0, 72.0, 87.0,
68.0, 80.0, 85.0, 90.0, 89.0,
66.0, 59.0, 60.0, 70.0, 71.0,
75.0, 76.0, 81.0, 85.0, 79.0,
90.0, 89.0, 97.0, 95.0, 100.0,
73.0, 76.0, 81.0, 79.0, 85.0,
44.0, 50.0, NA, 55.0, 53.0,
31.0, NA, 28.0, NA, NA,
60.0, 56.0, 67.0, 65.0, 72.0,
79.0, 68.0, 98.0, 88.0, 73.0,
92.0, 80.0, 91.0, NA, 78.0,
35.0, 62.0, 70.0, 74.0, 74.0,
39.0, 50.0, 71.0, NA, 75.0,
61.0, 88.0, 63.0, 83.0, 79.0,
70.0, 63.0, 64.0, 70.0, 91.0,
69.0, 53.0, 65.0, 66.0, NA,
60.0, 67.0, 75.0, 70.0, 80.0,
24.0, 48.0, 42.0, 52.0, 59.0,
40.0, 50.0, 66.0, 67.0, 59.0,
NA, 29.0, 49.0, 49.0, 45.0,
55.0, 63.0, 79.0, 69.0, 68.0,
60.0, 47.0, 52.0, 49.0, NA,
54.0, 64.0, 78.0, 68.0, 71.0,
61.0, 75.0, 77.0, 93.0, NA,
58.0, 47.0, 61.0, 50.0, 53.0,
40.0, 56.0, 51.0, 62.0, 57.0,
85.0, 74.0, 100.0, 92.0, 95.0,
60.0, 57.0, 83.0, 77.0, 82.0,
58.0, NA, 49.0, 55.0, 80.0,
60.0, 88.0, 89.0, 89.0, NA,
75.0, 73.0, 74.0, 82.0, 89.0,
54.0, 45.0, 43.0, 57.0, 50.0,
44.0, NA, 46.0, 73.0, 73.0,
30.0, 57.0, 30.0, 59.0, 40.0,
70.0, NA, 94.0, NA, 97.0,
51.0, 59.0, 49.0, 25.0, 61.0,
84.0, 89.0, 95.0, 93.0, 95.0,
58.0, 67.0, NA, 83.0, 83.0,
61.0, 60.0, 65.0, 89.0, 91.0,
56.0, 51.0, 61.0, 60.0, 55.0,
65.0, 68.0, 72.0, 66.0, 62.0,
66.0, 73.0, NA, 69.0, NA,
58.0, 58.0, 60.0, 74.0, 74.0,
39.0, 53.0, NA, NA, 65.0,
59.0, 76.0, 50.0, NA, 59.0,
60.0, NA, 60.0, NA, 57.0,
62.0, 74.0, 72.0, NA, 68.0,
NA, 66.0, 49.0, 66.0, NA,
85.0, 88.0, 91.0, 97.0, 96.0,
58.0, 49.0, 54.0, 39.0, 50.0,
50.0, 62.0, 50.0, 74.0, 64.0,
78.0, 93.0, NA, 93.0, 92.0,
58.0, 75.0, 61.0, 52.0, 68.0,
63.0, 59.0, 83.0, 67.0, 59.0,
67.0, 83.0, 90.0, 92.0, 93.0,
79.0, 78.0, 86.0, 94.0, 96.0,
49.0, 58.0, 66.0, 56.0, 68.0,
50.0, 70.0, 64.0, 68.0, 74.0,
48.0, 60.0, 56.0, 52.0, 55.0,
67.0, 74.0, NA, NA, NA,
95.0, 96.0, 95.0, 99.0, 98.0,
NA, 86.0, 88.0, 93.0, NA,
66.0, 71.0, 72.0, NA, NA,
74.0, NA, 88.0, 89.0, NA,
55.0, 42.0, 27.0, 41.0, 57.0,
63.0, 50.0, 52.0, 69.0, 67.0,
NA, 59.0, 75.0, 65.0, 62.0,
76.0, 83.0, 80.0, 86.0, 77.0,
67.0, 64.0, 63.0, 92.0, 72.0,
69.0, 70.0, 73.0, 76.0, 71.0,
56.0, NA, 43.0, 68.0, 80.0,
71.0, 55.0, 72.0, 66.0, 68.0,
60.0, 37.0, 88.0, 78.0, 85.0,
69.0, 66.0, 61.0, 73.0, 56.0,
71.0, 85.0, 74.0, 70.0, 96.0,
62.0, 65.0, NA, 22.0, 63.0,
53.0, 65.0, 67.0, 73.0, 71.0,
NA, 41.0, 34.0, NA, 49.0,
NA, 72.0, 66.0, 73.0, 75.0,
51.0, 63.0, 60.0, 67.0, NA,
56.0, 68.0, 75.0, NA, 90.0,
52.0, 62.0, 70.0, 60.0, 71.0,
16.0, 20.0, 23.0, 26.0, 33.0,
60.0, 54.0, 62.0, 77.0, 68.0,
88.0, 96.0, 89.0, NA, 98.0,
56.0, 65.0, 75.0, 83.0, 94.0,
85.0, 89.0, 92.0, NA, NA,
66.0, 70.0, 66.0, 69.0, 68.0,
72.0, 73.0, 74.0, 60.0, 54.0,
68.0, 82.0, 79.0, 77.0, 91.0,
60.0, 67.0, 66.0, 81.0, 82.0,
55.0, 47.0, 34.0, 65.0, 60.0,
80.0, 90.0, 89.0, 64.0, 88.0,
67.0, 87.0, 99.0, 87.0, 100.0,
79.0, 88.0, 73.0, 92.0, 90.0,
81.0, 78.0, 82.0, 87.0, 89.0,
65.0, 85.0, 75.0, 88.0, 73.0,
28.0, 2.0, 54.0, 61.0, 66.0,
88.0, 87.0, 90.0, 97.0, 96.0,
75.0, 76.0, 83.0, 87.0, 82.0,
65.0, 63.0, 69.0, 63.0, NA,
75.0, 73.0, 79.0, 81.0, 75.0,
70.0, 72.0, 85.0, 89.0, 93.0,
58.0, 46.0, 53.0, 70.0, 67.0,
83.0, 87.0, 86.0, 97.0, 87.0,
66.0, 69.0, 68.0, 89.0, 92.0,
65.0, 59.0, 66.0, 73.0, NA,
37.0, 36.0, NA, 43.0, 45.0,
78.0, 83.0, 63.0, 74.0, 88.0,
47.0, 55.0, 58.0, 67.0, 76.0,
97.0, 50.0, 45.0, 100.0, NA,
52.0, 56.0, 56.0, 100.0, 76.0,
44.0, 23.0, 39.0, 100.0, 56.0,
46.0, NA, 49.0, NA, 56.0,
55.0, 66.0, 66.0, 68.0, 70.0,
47.0, 52.0, 63.0, 66.0, 69.0,
49.0, 63.0, 59.0, 62.0, 62.0,
52.0, 50.0, 56.0, 62.0, 59.0,
50.0, 47.0, 54.0, 63.0, 63.0,
56.0, 69.0, 50.0, 75.0, 68.0,
73.0, 86.0, 75.0, 75.0, 82.0,
53.0, 57.0, 56.0, 69.0, 71.0,
70.0, 72.0, 70.0, 72.0, 72.0,
57.0, 55.0, 55.0, 53.0, 57.0,
51.0, 58.0, 64.0, 70.0, 63.0,
53.0, 57.0, 60.0, 70.0, 71.0,
50.0, 70.0, 70.0, 72.0, 81.0,
56.0, 55.0, 65.0, 70.0, 68.0,
49.0, 47.0, 89.0, 69.0, 69.0,
65.0, 72.0, 65.0, 76.0, 71.0,
75.0, 78.0, 92.0, 94.0, 94.0,
68.0, 71.0, 56.0, 75.0, 67.0,
55.0, 61.0, 60.0, 64.0, 62.0,
63.0, 71.0, 66.0, 67.0, 3.0,
56.0, 50.0, 0.0, 59.0, 82.0,
62.0, 67.0, 63.0, 57.0, 3.0,
62.0, 53.0, 73.0, 65.0, 75.0,
56.0, 53.0, 51.0, NA, 57.0,
92.0, 70.0, 71.0, 75.0, 71.0,
NA, 59.0, 65.0, 64.0, 49.0,
71.0, 69.0, 69.0, 75.0, 93.0,
86.0, 82.0, 90.0, 91.0, 94.0,
86.0, 87.0, 92.0, 94.0, 93.0,
NA, 36.0, 41.0, 46.0, NA,
58.0, 54.0, 65.0, 64.0, 62.0,
69.0, 71.0, 87.0, 81.0, 94.0,
84.0, 35.0, 29.0, 26.0, 96.0,
50.0, 70.0, 68.0, 80.0, 100.0,
76.0, 59.0, 73.0, 89.0, NA,
75.0, 90.0, 92.0, 97.0, 99.0,
43.0, 43.0, 71.0, 46.0, 44.0,
40.0, 44.0, 40.0, 40.0, 41.0,
48.0, 48.0, 80.0, 84.0, 82.0,
70.0, 73.0, 91.0, 71.0, NA,
46.0, 45.0, 49.0, 63.0, 65.0,
38.0, 58.0, 66.0, 65.0, 62.0,
49.0, 46.0, 66.0, 65.0, NA,
74.0, NA, 70.0, 85.0, 68.0,
46.0, 75.0, 78.0, 74.0, 68.0,
70.0, 71.0, 87.0, 91.0, 90.0,
65.0, 40.0, 60.0, 70.0, 72.0,
61.0, 54.0, 56.0, 55.0, 65.0,
48.0, 43.0, 44.0, 80.0, 61.0,
35.0, 55.0, 57.0, 34.0, 54.0,
73.0, 72.0, 84.0, 84.0, 87.0,
NA, 39.0, 62.0, 33.0, 38.0,
63.0, 93.0, 81.0, 89.0, 72.0,
73.0, 79.0, 71.0, 75.0, 92.0,
63.0, 84.0, 75.0, 64.0, 66.0,
NA, 66.0, 55.0, 66.0, 50.0,
66.0, 68.0, 69.0, 68.0, 71.0,
61.0, 57.0, 38.0, 43.0, 71.0,
63.0, 53.0, 90.0, 85.0, 87.0,
NA, 54.0, NA, 70.0, 51.0,
66.0, 70.0, 60.0, 69.0, 47.0,
73.0, 94.0, 92.0, 95.0, 96.0,
61.0, 50.0, 49.0, 50.0, NA,
91.0, NA, 94.0, NA, 99.0,
83.0, 76.0, 72.0, NA, 76.0,
51.0, 57.0, NA, NA, 60.0,
50.0, 60.0, 63.0, 50.0, 50.0,
64.0, 75.0, 88.0, NA, 91.0,
24.0, NA, NA, 45.0, NA,
NA, 91.0, 93.0, 96.0, 95.0,
51.0, 49.0, 49.0, 99.0, 72.0,
62.0, 89.0, 91.0, 90.0, 92.0,
56.0, 65.0, 69.0, 4.0, 71.0,
63.0, 60.0, 69.0, 64.0, 68.0,
49.0, 67.0, 77.0, 90.0, 70.0,
68.0, 71.0, 71.0, 79.0, 83.0,
59.0, 68.0, 64.0, 62.0, 63.0,
55.0, 59.0, 68.0, 71.0, 50.0,
80.0, 81.0, 76.0, 99.0, 88.0,
48.0, 46.0, 45.0, 66.0, 61.0,
60.0, 68.0, 69.0, 70.0, 73.0,
57.0, 64.0, 68.0, 68.0, 67.0,
63.0, 73.0, 70.0, 57.0, 63.0,
51.0, 60.0, 53.0, 51.0, 47.0,
60.0, 74.0, 69.0, 83.0, 69.0,
53.0, 58.0, 61.0, 61.0, 63.0,
60.0, 62.0, 57.0, 60.0, 73.0,
66.0, 75.0, 75.0, NA, 73.0,
53.0, 38.0, 36.0, NA, 40.0,
47.0, 37.0, 99.0, 98.0, 71.0,
60.0, 60.0, 94.0, NA, 62.0,
69.0, NA, 71.0, 87.0, 93.0,
48.0, 23.0, 26.0, 41.0, 52.0,
54.0, 38.0, 53.0, NA, 59.0,
55.0, 59.0, 57.0, 74.0, 76.0,
65.0, 80.0, 68.0, 72.0, 83.0,
38.0, 54.0, 44.0, 52.0, 54.0,
38.0, 28.0, 54.0, 66.0, 67.0,
17.0, 40.0, 8.0, 22.0, 25.0,
25.0, 26.0, 70.0, 18.0, 63.0,
NA, 73.0, 62.0, 85.0, 63.0,
48.0, 44.0, 59.0, 53.0, 61.0,
58.0, 70.0, 69.0, 61.0, 60.0,
80.0, 88.0, 87.0, 93.0, 93.0,
64.0, 66.0, 86.0, 86.0, 82.0,
23.0, 1.0, NA, 55.0, 52.0,
65.0, 68.0, 64.0, 69.0, NA,
85.0, 73.0, 58.0, 64.0, 60.0,
55.0, 67.0, 69.0, 77.0, 73.0,
45.0, 52.0, 56.0, 52.0, 75.0,
31.0, 19.0, 19.0, 50.0, 28.0,
20.0, 23.0, 38.0, 38.0, 50.0,
56.0, 68.0, 72.0, 70.0, 67.0,
64.0, 78.0, 80.0, 74.0, NA,
60.0, 54.0, 66.0, 67.0, NA,
44.0, 38.0, 47.0, 65.0, NA,
60.0, 61.0, 63.0, 60.0, 68.0,
54.0, 62.0, 69.0, 51.0, 65.0,
49.0, 47.0, 68.0, 62.0, 72.0,
44.0, 49.0, 47.0, 47.0, 52.0,
63.0, 73.0, 72.0, 88.0, 87.0,
65.0, 77.0, 93.0, NA, 71.0,
92.0, 93.0, 97.0, 97.0, 100.0,
60.0, 75.0, 83.0, 69.0, 68.0,
59.0, 68.0, 73.0, 75.0, 81.0,
83.0, 96.0, 88.0, 93.0, 96.0,
20.0, 61.0, 58.0, 21.0, 40.0,
76.0, 81.0, 85.0, 67.0, 69.0,
73.0, 49.0, 56.0, 42.0, 53.0,
NA, 46.0, 59.0, 57.0, NA,
69.0, 82.0, NA, 72.0, 71.0,
65.0, 69.0, NA, 79.0, 83.0,
42.0, 44.0, 47.0, 50.0, 59.0,
61.0, 68.0, 72.0, 73.0, 76.0,
87.0, 90.0, 92.0, 97.0, 90.0,
67.0, 74.0, 69.0, 77.0, 95.0,
90.0, 92.0, 90.0, 91.0, 100.0,
65.0, 59.0, 51.0, NA, 72.0,
NA, 59.0, 61.0, 100.0, 62.0,
85.0, 90.0, 89.0, 90.0, 88.0,
70.0, 69.0, 64.0, 78.0, 87.0,
59.0, 75.0, 69.0, 70.0, 76.0,
71.0, 65.0, 91.0, 94.0, NA,
76.0, 80.0, 73.0, 75.0, 78.0,
69.0, NA, 80.0, 71.0, 96.0,
48.0, 57.0, 60.0, 59.0, 67.0,
54.0, 41.0, 59.0, 68.0, 38.0,
58.0, 67.0, 94.0, 95.0, 93.0,
54.0, 48.0, 58.0, 70.0, 70.0,
64.0, 54.0, 52.0, 35.0, 51.0,
53.0, 72.0, 68.0, 66.0, 79.0,
60.0, 70.0, 67.0, 65.0, 72.0,
31.0, NA, 40.0, 46.0, 52.0,
59.0, 61.0, 69.0, NA, NA,
59.0, 50.0, 100.0, 61.0, 97.0,
67.0, 76.0, 74.0, 75.0, 83.0,
62.0, 92.0, NA, NA, 55.0,
52.0, 66.0, 92.0, 76.0, 64.0,
42.0, 6.0, 48.0, 75.0, 1.0,
81.0, 73.0, NA, NA, 72.0,
34.0, 36.0, 44.0, 54.0, 40.0,
66.0, 70.0, NA, NA, 91.0,
86.0, 56.0, 70.0, NA, 46.0,
66.0, 72.0, 79.0, NA, 79.0,
49.0, 73.0, 70.0, 78.0, 68.0,
62.0, 50.0, 61.0, 81.0, 65.0,
93.0, 90.0, 97.0, 98.0, 100.0,
22.0, 20.0, NA, 27.0, 68.0,
70.0, 70.0, 79.0, 83.0, 85.0,
81.0, NA, NA, 91.0, 92.0,
58.0, 73.0, 73.0, 97.0, 73.0,
NA, 5.0, 32.0, 39.0, 44.0,
95.0, 98.0, 97.0, 97.0, 95.0,
54.0, 60.0, 60.0, 70.0, 71.0,
72.0, 70.0, NA, NA, 72.0,
NA, 49.0, 61.0, 56.0, 54.0,
NA, 82.0, 72.0, 95.0, 65.0,
NA, 71.0, 87.0, 67.0, 76.0,
68.0, 56.0, 88.0, 80.0, NA,
55.0, 88.0, 76.0, 94.0, 89.0,
78.0, 77.0, 88.0, 88.0, 87.0,
86.0, 91.0, 91.0, 88.0, 96.0,
62.0, 69.0, NA, 49.0, NA,
31.0, 45.0, 50.0, 41.0, 58.0,
46.0, 64.0, NA, 73.0, 70.0,
61.0, 67.0, 64.0, 65.0, 58.0,
59.0, 56.0, 59.0, 75.0, 56.0,
NA, 66.0, 57.0, 71.0, NA,
29.0, 44.0, 44.0, 42.0, 63.0,
49.0, 52.0, 54.0, 62.0, 59.0,
43.0, 78.0, 70.0, 90.0, 69.0,
83.0, 68.0, 76.0, 71.0, 69.0,
63.0, 65.0, 79.0, 72.0, 65.0,
55.0, 59.0, 59.0, 58.0, 59.0,
51.0, 69.0, 52.0, 60.0, 38.0,
28.0, 61.0, NA, 59.0, NA,
58.0, 56.0, 67.0, 66.0, 67.0,
54.0, 47.0, 58.0, 59.0, 58.0,
NA, 40.0, 55.0, 68.0, 70.0,
46.0, 61.0, 54.0, 70.0, 73.0,
44.0, 56.0, 66.0, 54.0, 67.0,
63.0, 69.0, 66.0, 68.0, 71.0,
46.0, 54.0, 64.0, 49.0, 68.0,
NA, 69.0, 78.0, 58.0, 64.0,
NA, 57.0, 60.0, 74.0, 51.0,
68.0, 37.0, 51.0, 56.0, 56.0,
74.0, 64.0, 78.0, 80.0, 73.0,
70.0, 65.0, 65.0, 59.0, 66.0,
83.0, 73.0, 91.0, 93.0, 92.0,
74.0, 70.0, 53.0, 64.0, 64.0,
59.0, 56.0, 80.0, 57.0, 54.0,
72.0, NA, NA, 49.0, 50.0,
78.0, 91.0, 89.0, 87.0, 82.0,
31.0, 1.0, 47.0, 49.0, 57.0,
84.0, 77.0, 65.0, 87.0, 64.0,
63.0, 0.0, 64.0, 71.0, 73.0,
42.0, 50.0, 57.0, 57.0, 62.0,
31.0, 37.0, 57.0, 41.0, 39.0,
54.0, 71.0, 88.0, 88.0, 87.0,
54.0, NA, 53.0, 46.0, 71.0,
60.0, 62.0, 62.0, 61.0, 87.0,
81.0, 61.0, 72.0, 80.0, 88.0,
88.0, 90.0, 70.0, 92.0, 93.0,
24.0, 30.0, 27.0, 27.0, 20.0,
63.0, 58.0, 64.0, 51.0, 64.0,
79.0, 80.0, 82.0, 87.0, 86.0,
46.0, 48.0, 63.0, 68.0, 65.0,
66.0, 69.0, 76.0, 69.0, 65.0,
68.0, 69.0, 68.0, 70.0, 87.0,
59.0, 60.0, 64.0, 61.0, 72.0,
73.0, 86.0, 88.0, 90.0, NA,
69.0, 65.0, 69.0, 74.0, NA),byrow=T, ncol=5)
##
library(moments)
apply(y, 2, mean, na.rm=T)
apply(y, 2, sd, na.rm=T)
apply(y, 2, skewness, na.rm=T)
apply(y, 2, kurtosis, na.rm=T)
##
par(mfrow=c(3,2))
qqnorm(y[,1], main='')
qqline(y[,1])
qqnorm(y[,2], main='')
qqline(y[,2])
qqnorm(y[,3], main='')
qqline(y[,3])
qqnorm(y[,4], main='')
qqline(y[,4])
qqnorm(y[,5], main='')
qqline(y[,5])
## boxplots
par(family='times')
boxplot(y[,1],y[,2],y[,3],y[,4],y[,5], names=c('7','8','9','10','11'), ylab='Math Score', xlab='Grade')
boxplot(rnorm(375, mean(y[,1],na.rm=T), sd(y[,1],na.rm=T)), rnorm(375, mean(y[,2],na.rm=T), sd(y[,2],na.rm=T)),rnorm(375, mean(y[,3],na.rm=T), sd(y[,3],na.rm=T)),rnorm(375, mean(y[,4],na.rm=T), sd(y[,4],na.rm=T)),rnorm(375, mean(y[,5],na.rm=T), sd(y[,5],na.rm=T)),add=T, col='red')