5 Matrices and efficiency measures for block designs

This function returns the point-block incidence matrix N of the block design D. This matrix has rows indexed by the points of D and columns by the blocks of D, with the (i,j)-entry of N being the number of times point i occurs in D.blocks[j].

The returned matrix N is immutable.

gap> D:=DualBlockDesign(AGPointFlatBlockDesign(2,3,1));;
gap> BlockDesignBlocks(D);
[ [ 1, 2, 3, 4 ], [ 1, 5, 6, 7 ], [ 1, 8, 9, 10 ], [ 2, 5, 8, 11 ], 
  [ 2, 7, 9, 12 ], [ 3, 5, 10, 12 ], [ 3, 6, 9, 11 ], [ 4, 6, 8, 12 ], 
  [ 4, 7, 10, 11 ] ]
gap> PointBlockIncidenceMatrix(D);
[ [ 1, 1, 1, 0, 0, 0, 0, 0, 0 ], [ 1, 0, 0, 1, 1, 0, 0, 0, 0 ], 
  [ 1, 0, 0, 0, 0, 1, 1, 0, 0 ], [ 1, 0, 0, 0, 0, 0, 0, 1, 1 ], 
  [ 0, 1, 0, 1, 0, 1, 0, 0, 0 ], [ 0, 1, 0, 0, 0, 0, 1, 1, 0 ], 
  [ 0, 1, 0, 0, 1, 0, 0, 0, 1 ], [ 0, 0, 1, 1, 0, 0, 0, 1, 0 ], 
  [ 0, 0, 1, 0, 1, 0, 1, 0, 0 ], [ 0, 0, 1, 0, 0, 1, 0, 0, 1 ], 
  [ 0, 0, 0, 1, 0, 0, 1, 0, 1 ], [ 0, 0, 0, 0, 1, 1, 0, 1, 0 ] ]

This function returns the concurrence matrix L of the block design D. This matrix is equal to NN^T, where N is the point-block incidence matrix of D (see PointBlockIncidenceMatrix) and N^T is the transpose of N. If D is a binary block design then the (i,j)-entry of its concurrence matrix is the number of blocks containing points i and j.

The returned matrix L is immutable.

gap> D:=DualBlockDesign(AGPointFlatBlockDesign(2,3,1));;
gap> BlockDesignBlocks(D);
[ [ 1, 2, 3, 4 ], [ 1, 5, 6, 7 ], [ 1, 8, 9, 10 ], [ 2, 5, 8, 11 ], 
  [ 2, 7, 9, 12 ], [ 3, 5, 10, 12 ], [ 3, 6, 9, 11 ], [ 4, 6, 8, 12 ], 
  [ 4, 7, 10, 11 ] ]
gap> ConcurrenceMatrix(D);
[ [ 3, 1, 1, 1, 1, 1, 1, 1, 1, 1, 0, 0 ], 
  [ 1, 3, 1, 1, 1, 0, 1, 1, 1, 0, 1, 1 ], 
  [ 1, 1, 3, 1, 1, 1, 0, 0, 1, 1, 1, 1 ], 
  [ 1, 1, 1, 3, 0, 1, 1, 1, 0, 1, 1, 1 ], 
  [ 1, 1, 1, 0, 3, 1, 1, 1, 0, 1, 1, 1 ], 
  [ 1, 0, 1, 1, 1, 3, 1, 1, 1, 0, 1, 1 ], 
  [ 1, 1, 0, 1, 1, 1, 3, 0, 1, 1, 1, 1 ], 
  [ 1, 1, 0, 1, 1, 1, 0, 3, 1, 1, 1, 1 ], 
  [ 1, 1, 1, 0, 0, 1, 1, 1, 3, 1, 1, 1 ], 
  [ 1, 0, 1, 1, 1, 0, 1, 1, 1, 3, 1, 1 ], 
  [ 0, 1, 1, 1, 1, 1, 1, 1, 1, 1, 3, 0 ], 
  [ 0, 1, 1, 1, 1, 1, 1, 1, 1, 1, 0, 3 ] ]

This function returns the information matrix C of the block design D.

This matrix is defined as follows. Suppose D has v points and b blocks, let R be the vtimesv diagonal matrix whose (i,i)-entry is the replication number of the point i, let N be the point-block incidence matrix of D (see PointBlockIncidenceMatrix), and let K be the btimesb diagonal matrix whose (j,j)-entry is the length of D.blocks[j]. Then the information matrix of D is C:=R-NK^-1N^T. If D is a 1-(v,k,r) design then this expression for C simplifies to rI-k^-1L, where I is the vtimesv identity matrix and L is the concurrence matrix of D (see ConcurrenceMatrix).

The returned matrix C is immutable.

gap> D:=DualBlockDesign(AGPointFlatBlockDesign(2,3,1));;
gap> BlockDesignBlocks(D);
[ [ 1, 2, 3, 4 ], [ 1, 5, 6, 7 ], [ 1, 8, 9, 10 ], [ 2, 5, 8, 11 ], 
  [ 2, 7, 9, 12 ], [ 3, 5, 10, 12 ], [ 3, 6, 9, 11 ], [ 4, 6, 8, 12 ], 
  [ 4, 7, 10, 11 ] ]
gap> InformationMatrix(D);
[ [ 9/4, -1/4, -1/4, -1/4, -1/4, -1/4, -1/4, -1/4, -1/4, -1/4, 0, 0 ], 
  [ -1/4, 9/4, -1/4, -1/4, -1/4, 0, -1/4, -1/4, -1/4, 0, -1/4, -1/4 ], 
  [ -1/4, -1/4, 9/4, -1/4, -1/4, -1/4, 0, 0, -1/4, -1/4, -1/4, -1/4 ], 
  [ -1/4, -1/4, -1/4, 9/4, 0, -1/4, -1/4, -1/4, 0, -1/4, -1/4, -1/4 ], 
  [ -1/4, -1/4, -1/4, 0, 9/4, -1/4, -1/4, -1/4, 0, -1/4, -1/4, -1/4 ], 
  [ -1/4, 0, -1/4, -1/4, -1/4, 9/4, -1/4, -1/4, -1/4, 0, -1/4, -1/4 ], 
  [ -1/4, -1/4, 0, -1/4, -1/4, -1/4, 9/4, 0, -1/4, -1/4, -1/4, -1/4 ], 
  [ -1/4, -1/4, 0, -1/4, -1/4, -1/4, 0, 9/4, -1/4, -1/4, -1/4, -1/4 ], 
  [ -1/4, -1/4, -1/4, 0, 0, -1/4, -1/4, -1/4, 9/4, -1/4, -1/4, -1/4 ], 
  [ -1/4, 0, -1/4, -1/4, -1/4, 0, -1/4, -1/4, -1/4, 9/4, -1/4, -1/4 ], 
  [ 0, -1/4, -1/4, -1/4, -1/4, -1/4, -1/4, -1/4, -1/4, -1/4, 9/4, 0 ], 
  [ 0, -1/4, -1/4, -1/4, -1/4, -1/4, -1/4, -1/4, -1/4, -1/4, 0, 9/4 ] ]

5.2 The function BlockDesignEfficiency

Let D be a 1-(v,k,r) design with v>1, let eps be a positive rational number (default: 10^-6), and let includeMV be a boolean (default: false). Then this function returns a record eff containing information on statistical efficiency measures of D. These measures are defined below. See Extrep, BaCa and BaRo for further details. All returned results are computed using exact algebraic computation.

The component eff.A contains the A-efficiency measure for D, eff.Dpowered contains the D-efficiency measure of D raised to the power v-1, and eff.Einterval is a list [a,b] of non-negative rational numbers such that if x is the E-efficiency measure of D then alexleb, b-aleeps, and if x is rational then a=x=b. Moreover eff.CEFpolynomial contains the monic polynomial over the rationals whose zeros (counting multiplicities) are the canonical efficiency factors of the design D. If includeMV=true then additional work is done to compute the MV- (also called E'-) efficiency measure, and then eff.MV contains the value of this measure. (This component may be set even if includeMV=false, as a byproduct of other computation.)

We now define the canonical efficiency factors and the A-, D-, E-, and MV-efficiency measures of a 1-design.

Let D be a 1-(v,k,r) design with vge2, let C be the information matrix of D (see InformationMatrix), and let F:=r^-1C. The eigenvalues of F are all real and lie in the interval [0,1]. At least one of these eigenvalues is zero: an associated eigenvector is the all-1 vector. The remaining eigenvalues delta₁ledelta₂le cdotsledelta_v-1 of F are called the canonical efficiency factors of D. These are all non-zero if and only if D is connected (that is, the point-block incidence graph of D is a connected graph).

If D is not connected, then the A-, D-, E-, and MV-efficiency measures of D are all defined to be zero. Otherwise, the A-efficiency measure is (v-1)/sum_i=1^v-11/delta_i (the harmonic mean of the canonical efficiency factors), the D-efficiency measure is (prod_i=1^v-1delta_i)^1/(v-1) (the geometric mean of the canonical efficiency factors), and the E-efficiency measure is delta₁ (the minimum of the canonical efficiency factors).

If D is connected, and the MV-efficiency measure is required, then it is computed as follows. Let F:=r^-1C be as before, and let P:=v^-1J, where J is the vtimesv all-1 matrix. Set M:=(F+P)^-1-P, making M the ``Moore-Penrose inverse'' of F (see BaCa). Then the MV-efficiency measure of D is the minimum value (over all i,jin{1,...,v}, inot=j) of 2/(M_ii+M_jj-M_ij-M_ji).

gap> D:=DualBlockDesign(AGPointFlatBlockDesign(2,3,1));;
gap> BlockDesignBlocks(D);
[ [ 1, 2, 3, 4 ], [ 1, 5, 6, 7 ], [ 1, 8, 9, 10 ], [ 2, 5, 8, 11 ], 
  [ 2, 7, 9, 12 ], [ 3, 5, 10, 12 ], [ 3, 6, 9, 11 ], [ 4, 6, 8, 12 ], 
  [ 4, 7, 10, 11 ] ]
gap> BlockDesignEfficiency(D);
rec( A := 33/41, 
  CEFpolynomial := x_1^11-9*x_1^10+147/4*x_1^9-719/8*x_1^8+18723/128*x_1^7-106\
47/64*x_1^6+138159/1024*x_1^5-159813/2048*x_1^4+2067201/65536*x_1^3-556227/655\
36*x_1^2+89667/65536*x_1-6561/65536, Dpowered := 6561/65536, 
  Einterval := [ 3/4, 3/4 ] )
gap> BlockDesignEfficiency(D,10^(-4),true);
rec( A := 33/41, 
  CEFpolynomial := x_1^11-9*x_1^10+147/4*x_1^9-719/8*x_1^8+18723/128*x_1^7-106\
47/64*x_1^6+138159/1024*x_1^5-159813/2048*x_1^4+2067201/65536*x_1^3-556227/655\
36*x_1^2+89667/65536*x_1-6561/65536, Dpowered := 6561/65536, 
  Einterval := [ 3/4, 3/4 ], MV := 3/4 )

5.3 Computing an interval for a certain real zero of a rational polynomial

We document a DESIGN package utility function used in the calculation of the Einterval component above, but is more widely applicable.

Suppose that f is a univariate polynomial over the rationals, a, b are rational numbers with aleb, and eps is a positive rational number.

If f has no real zero in the closed interval [a,b], then this function returns the empty list. Otherwise, let alpha be the least real zero of f such that alealphaleb. Then this function returns a list [c,d] of rational numbers, with clealphaled and d-cleeps. Moreover, c=d if and only if alpha is rational (in which case alpha=c=d).

gap> x:=Indeterminate(Rationals,1);                     
x_1
gap> f:=(x+3)*(x^2-3);
x_1^3+3*x_1^2-3*x_1-9
gap> L:=DESIGN_IntervalForLeastRealZero(f,-5,5,10^(-3));
[ -3, -3 ]
gap> L:=DESIGN_IntervalForLeastRealZero(f,-2,5,10^(-3));
[ -14193/8192, -7093/4096 ]
gap> List(L,Float);             
[ -1.73254, -1.73169 ]
gap> L:=DESIGN_IntervalForLeastRealZero(f,0,5,10^(-3));
[ 14185/8192, 7095/4096 ]
gap> List(L,Float);           
[ 1.73157, 1.73218 ]
gap> L:=DESIGN_IntervalForLeastRealZero(f,0,5,10^(-5));
[ 454045/262144, 908095/524288 ]
gap> List(L,Float);                  
[ 1.73204, 1.73205 ]
gap> L:=DESIGN_IntervalForLeastRealZero(f,2,5,10^(-5));
[  ]

[Up] [Previous] [Next] [Index]