# 编程代写 Math 558 Lecture #29 – cscodehelp代写

Math 558 Lecture #29

More on BIBDs in R The package “ibd”
In this lecture we will use the package “ibd”.

Copyright By cscodehelp代写 加微信 cscodehelp

help(ibd) will provide a lot of information about the available functions. I will discuss some in today’s lecture. To start with we will install the package first 1
install.packages(“ibd”)
library(ibd)
Design1 <- bibd(v = 7, b = 7, r = 3, k = 3, lambda = 1) Here the symbol v is used for the treatments. Hence we can obtain the plan for any BIBD if it exists for the given parameters2. 1it was not easy for me to install. I had to update some packages and download and install XQuartz-2.8.1.dmg 2a BIBD with given parameters may not exist. We will discuss more about that in the coming lectures More on BIBDs in R The design that we will obtain is given below. [,1] [,2] [,3] Block-1 1 Block-2 5 Block-3 3 Block-4 1 Block-5 1 Block-6 2 Block-7 2 Concurrence matrix for the generated design \$NNP(Concurrence matrix) [,1] [,2] [,3] [,4] [,5] [,6] [,7] [1,] 3 1 1 1 1 1 1 [2,] 1 3 1 1 1 1 1 [3,] 1 1 3 1 1 1 1 [4,] 1 1 1 3 1 1 1 [5,] 1 1 1 1 3 1 1 [6,] 1 1 1 1 1 3 1 [7,] 1 1 1 1 1 1 3 Concurrence matrix of a design The matrix given in the previous slide is called the concurrence matrix of the design. This matrix gives the total number of occurrences of all the treatments in the diagonal and pair wise concurrences in the off diagonal places. We can denote this matrix by M. The general form of the matrix is given below This a t × t matrix r λ λ .... λ λ r λ .... λ ..... ...... ..... .... .... M = ..... ...... ..... .... ...  ..... ...... .... λ r  Incidence matrix Another matrix that you will see in the R output is the incidence matrix given below 1 0 0 1 1 0 0 0 0 0 0 1 1 1 0 0 1 1 0 0 1 N=1 0 1 0 0 1 0 1 1 0 0 0 0 1 0 1 0 1 0 1 0 Incidence matrix The incident matrix is given by N = (nij)t×b 􏱁1 if treatment i is in block j 0 otherwise It can be seen that M = NNT Note: It is the design in which a treatment occurs at most once in a block. A BIBD is always a binary design. yij = μ + bi + τj + εij which is identical to the model for a randomized complete block design given . However, the analysis is slightly different due to the missing observations as the blocks are incomplete. Hereyij is the ith observation in the jth block, μi is the overall mean, τi is the effect of the ith treatment, bj is the effect of the jth block, and εij random error component. Also, εij ∼iid N(0, σ2) Statistical Analysis The total variability in the data is given by SST 2 y2.. SST = ∑∑yij − N ij The total variability can be partition as SST = SStr(adjusted) + SSB + SSE The sum of squares for treatments are adjusted each treatment appears in a different set of r blocks. Therefore the difference between the treatment totals (unadjusted), yi. for i = 1, 2, ..t are effected by the difference in blocks as well. The block sum of squares are 1k y2 SSB= ∑y2.j−.. Statistical Analysis y.j = y1j + y2j + ..ykj is the total in jth block. Sum of squares for blocks has b − 1 degrees of freedom. Define the quantity Qi as Qi =yi.−k∑nijy.j i=1,2..t Where nij = 1 if treatment i appears in block j and 0 otherwise The adjusted treatment sum of squares is SSTr(adjusted) = k ∑ti=1 Q2i It has t − 1 degrees of freedom. The SSE has N − 1 − (a − 1) − (b − 1) = N − a − b + 1 degrees of freedom. The appropriate test statistic is F = MSTr(adjusted) MSE Balanced Incomplete Block Design for Catalyst Experiment Treatment Block (Batch of Raw Material) 1 2 3 4 yi. (Catalyst) 1 73 74 — 71 218 2 — 75 67 72 214 3 73 75 68 — 216 4 75 — 72 75 222 y.j 221 224 207 218 870 y. Statistical Analysis: Catalyst example This is a BIBD with t = 4, b = 4, k = r = 3, λ = 2, and N = 12 experimental units. 2 y2.. 8702 SST=∑∑yij−N=63.156− 12 =81.00 1k y2 SSB= ∑y2.j− .. =1[(221)2 + (207)2 + (224)2 + (218)2] − 8702 3 12 =55.00 Statistical Analysis Q1 = (218) − 1/3(221 + 224 + 218) = −9/3 Q2 =−7/3 Q3 =−4/3 Q4 =20/3 SSTr(adjusted) = k ∑ti=1 Q2i λt = 3[(−9/3)2 + (−7/3)2 + (−4/3)2 + (20/3)2 ] 2×4 SSE = 81.00 − 22.75 − 55.00 = 3.25 F = 11.6, p − value = 0.107 ANOVA table Statistical analysis in R data = taste different codes Model1<-aov(score panelist + recipe, data = taste) > summary(Model1)
Df Sum Sq Mean Sq F value Pr(>F)
panelist 11 19.333 1.7576 2.301 0.1106
recipe 3 9.125 3.0417 3.982 0.0465 * Residuals 9 6.875 0.7639

Further ivestigation install.packages(“multcomp”) library(multcomp)
contr <- glht(Model1, linfct = mcp(recipe = "Tukey")) summary(contr, test = adjusted("none")) Estimate B-A==0 0.750 C-A==0 1.375 D-A==0 -0.625 C-B==0 0.625 D-B==0 -1.375 D-C==0 -2.000 Std. Error 0.618 0.618 0.618 0.618 0.618 0.618 t value 1.214 2.225 -1.011 1.011 -2.225 -3.236 Pr(>|t|) 0.2558 0.0531 0.3383 0.3383 0.0531 . 0.0102 *