In [1]:

    

alpha <- 0.05


    

In [2]:

    

data <- data.frame(B1=c(10, 11, 6), B2=c(21, 27, 19), B3=c(15, 21, 27), B4=c(6, 13, 24),
                   row.names=c("A1", "A2", "A3"))
data


    

    
Out[2]:





B1
B2
B3
B4

	
A1
10
21
15
6
	
A2
11
27
21
13
	
A3
6
19
27
24

In [3]:

    

a <- nrow(data)
b <- ncol(data)
n <- sum(data)
c(a, b, n)


    

    
Out[3]:





	
3
	
4
	
200

In [4]:

    

p_i <- 1/n * rowSums(data)
p_i


    

    
Out[4]:





	
A1
		
0.26
	
A2
		
0.36
	
A3
		
0.38

In [5]:

    

p_j <- 1/n * colSums(data)
p_j


    

    
Out[5]:





	
B1
		
0.135
	
B2
		
0.335
	
B3
		
0.315
	
B4
		
0.215

In [6]:

    

p_i <- 1/n * rowSums(data)
p_i


    

    
Out[6]:





	
A1
		
0.26
	
A2
		
0.36
	
A3
		
0.38

In [7]:

    

p_ij <- p_i %*% t(p_j)
p_ij


    

    
Out[7]:





B1
B2
B3
B4

	
0.0351
0.0871
0.0819
0.0559
	
0.0486
0.1206
0.1134
0.0774
	
0.0513
0.1273
0.1197
0.0817

In [8]:

    

Q_ij <- (data - n*p_ij)^2/(n * p_ij)
Q_ij


    

    
Out[8]:





B1
B2
B3
B4

	
A1
1.2650142
0.7357290
0.1162637
2.4000358
	
A2
0.1685597
0.3438806
0.1244444
0.3973127
	
A3
1.7687719
1.6391045
0.3911278
3.5909180

In [9]:

    

sum(Q_ij)


    

    
Out[9]:




12.941162399907

In [10]:

    

qchisq(1-alpha, (a-1)*(b-1))


    

    
Out[10]:




12.591587243744