Index Matrices

The indMatrix class is the class of row and column index matrices, stored as 1-based integer index vectors. A row (column) index matrix is a matrix whose rows (columns) are standard unit vectors. Such matrices are useful when mapping observations to discrete sets of covariate values.

Multiplying a matrix on the left by a row index matrix is equivalent to indexing its rows, i.e., sampling the rows “with replacement”. Analogously, multiplying a matrix on the right by a column index matrix is equivalent to indexing its columns. Indeed, such products are implemented in Matrix as indexing operations; see ‘Details’ below.

A matrix whose rows and columns are standard unit vectors is called a permutation matrix. This special case is designated by the pMatrix class, a direct subclass of indMatrix.

Objects from the Class

Objects can be created explicitly with calls of the form new("indMatrix", ...), but they are more commonly created by coercing 1-based integer index vectors, with calls of the form as(., "indMatrix"); see ‘Methods’ below.

Slots

margin: an integer, either 1 or 2, specifying whether the matrix is a row (1) or column (2) index.
perm: a 1-based integer index vector, i.e., a vector of length Dim[margin] with elements taken from 1:Dim[1+margin%%2].
Dim,Dimnames: inherited from virtual superclass Matrix.

Extends

Classes "sparseMatrix" and "generalMatrix", directly.

Methods

%*%: signature(x = "indMatrix", y = "Matrix") and others listed by showMethods("%*%", classes = "indMatrix"): matrix products implemented where appropriate as indexing operations.
coerce: signature(from = "numeric", to = "indMatrix"): supporting typical indMatrix construction from a vector of positive integers. Row indexing is assumed.
coerce: signature(from = "list", to = "indMatrix"): supporting indMatrix construction for row and column indexing, including index vectors of length 0 and index vectors whose maximum is less than the number of rows or columns being indexed.
coerce: signature(from = "indMatrix", to = "matrix"): coercion to a traditional matrix of logical type, with FALSE and TRUE in place of 0 and 1.
t: signature(x = "indMatrix"): the transpose, which is an indMatrix with identical perm but opposite margin.
rowSums,rowMeans,colSums,colMeans: signature(x = "indMatrix"): row and column sums and means.
rbind2,cbind2: signature(x = "indMatrix", y = "indMatrix"): row-wise catenation of two row index matrices with equal numbers of columns and column-wise catenation of two column index matrices with equal numbers of rows.
kronecker: signature(X = "indMatrix", Y = "indMatrix"): Kronecker product of two row index matrices or two column index matrices, giving the row or column index matrix corresponding to their “interaction”.

Author

Fabian Scheipl at uni-muenchen.de, building on the existing class pMatrix after a nice hike's conversation with Martin Maechler. Methods for crossprod(x, y) and kronecker(x, y) with both arguments inheriting from indMatrix were made considerably faster thanks to a suggestion by Boris Vaillant. Diverse tweaks by Martin Maechler and Mikael Jagan, notably the latter's implementation of margin, prior to which the indMatrix class was designated only for row index matrices.

Details

The transpose of an index matrix is an index matrix with identical perm but opposite margin. Hence the transpose of a row index matrix is a column index matrix, and vice versa.

The cross product of a row index matrix R and itself is a diagonal matrix whose diagonal entries are the the number of entries in each column of R.

Given a row index matrix R with perm slot p, a column index matrix C with perm slot q, and a matrix M with conformable dimensions, we have

\(R M\)	=	`R %*% M`	=	`M[p, ]`
\(M C\)	=	`M %*% C`	=	`M[, q]`
\(C'M\)	=	`crossprod(C, M)`	=	`M[q, ]`
\(MR'\)	=	`tcrossprod(M, R)`	=	`M[, p]`
\(R'R\)	=	`crossprod(R)`	=	`Diagonal(x=tabulate(p, ncol(R)))`
\(CC'\)	=	`tcrossprod(C)`	=	`Diagonal(x=tabulate(q, nrow(C)))`

Operations on index matrices that result in index matrices will accordingly return an indMatrix. These include products of two column index matrices and (equivalently) column-indexing of a column index matrix (when dimensions are not dropped). Most other operations on indMatrix treat them as sparse nonzero pattern matrices (i.e., inheriting from virtual class nsparseMatrix). Hence vector-valued subsets of indMatrix, such as those given by diag, are always of type "logical".

Examples

p1 <- as(c(2,3,1), "pMatrix")
(sm1 <- as(rep(c(2,3,1), e=3), "indMatrix"))
#> 9 x 3 sparse Matrix of class "indMatrix"
#>            
#>  [1,] . | .
#>  [2,] . | .
#>  [3,] . | .
#>  [4,] . . |
#>  [5,] . . |
#>  [6,] . . |
#>  [7,] | . .
#>  [8,] | . .
#>  [9,] | . .
stopifnot(all(sm1 == p1[rep(1:3, each=3),]))

## row-indexing of a <pMatrix> turns it into an <indMatrix>:
class(p1[rep(1:3, each=3),])
#> [1] "indMatrix"
#> attr(,"package")
#> [1] "Matrix"

set.seed(12) # so we know '10' is in sample
## random index matrix for 30 observations and 10 unique values:
(s10 <- as(sample(10, 30, replace=TRUE),"indMatrix"))
#> 30 x 10 sparse Matrix of class "indMatrix"
#>                          
#>  [1,] . | . . . . . . . .
#>  [2,] . . . . . . . . . |
#>  [3,] . . . . . . | . . .
#>  [4,] . . . . | . . . . .
#>  [5,] . . . . | . . . . .
#>  [6,] . | . . . . . . . .
#>  [7,] . . . . . . . | . .
#>  [8,] . | . . . . . . . .
#>  [9,] . | . . . . . . . .
#> [10,] . . . . . . . . | .
#> [11,] . . . . . | . . . .
#> [12,] . . . . . . . . . |
#> [13,] . . . . . . . . . |
#> [14,] . . . . . . | . . .
#> [15,] . . . . . . . | . .
#> [16,] . . . | . . . . . .
#> [17,] . . . . . . . | . .
#> [18,] . . . . . . . . | .
#> [19,] . . . . . . . | . .
#> [20,] . . . . . . . . . |
#> [21,] . . . . . . | . . .
#> [22,] . . . . | . . . . .
#> [23,] . . . . . . . . | .
#> [24,] . | . . . . . . . .
#> [25,] . . . . . . | . . .
#> [26,] . . . . . . | . . .
#> [27,] . . . . . . . . . |
#> [28,] . . . | . . . . . .
#> [29,] . . . | . . . . . .
#> [30,] . . | . . . . . . .

## Sample rows of a numeric matrix :
(mm <- matrix(1:10, nrow=10, ncol=3))
#>       [,1] [,2] [,3]
#>  [1,]    1    1    1
#>  [2,]    2    2    2
#>  [3,]    3    3    3
#>  [4,]    4    4    4
#>  [5,]    5    5    5
#>  [6,]    6    6    6
#>  [7,]    7    7    7
#>  [8,]    8    8    8
#>  [9,]    9    9    9
#> [10,]   10   10   10
s10 %*% mm
#> 30 x 3 Matrix of class "dgeMatrix"
#>       [,1] [,2] [,3]
#>  [1,]    2    2    2
#>  [2,]   10   10   10
#>  [3,]    7    7    7
#>  [4,]    5    5    5
#>  [5,]    5    5    5
#>  [6,]    2    2    2
#>  [7,]    8    8    8
#>  [8,]    2    2    2
#>  [9,]    2    2    2
#> [10,]    9    9    9
#> [11,]    6    6    6
#> [12,]   10   10   10
#> [13,]   10   10   10
#> [14,]    7    7    7
#> [15,]    8    8    8
#> [16,]    4    4    4
#> [17,]    8    8    8
#> [18,]    9    9    9
#> [19,]    8    8    8
#> [20,]   10   10   10
#> [21,]    7    7    7
#> [22,]    5    5    5
#> [23,]    9    9    9
#> [24,]    2    2    2
#> [25,]    7    7    7
#> [26,]    7    7    7
#> [27,]   10   10   10
#> [28,]    4    4    4
#> [29,]    4    4    4
#> [30,]    3    3    3

set.seed(27)
IM1 <- as(sample(1:20, 100, replace=TRUE), "indMatrix")
IM2 <- as(sample(1:18, 100, replace=TRUE), "indMatrix")
(c12 <- crossprod(IM1,IM2))
#> 20 x 18 sparse Matrix of class "dgTMatrix"
#>                                          
#>  [1,] . . 2 . . . . 1 1 . . . 1 . . . . 1
#>  [2,] . . . . . . 1 . 1 . 1 . . . . . 1 .
#>  [3,] . 1 . 1 1 1 . . . . . . . . 1 1 . 1
#>  [4,] . . . 1 1 . . . 1 . 1 . . 1 1 1 1 .
#>  [5,] . 1 1 . 1 . 1 . . . 1 . . . 1 2 . 1
#>  [6,] . . . 1 . . . . . 1 . . . 1 . . . .
#>  [7,] . . . . . . . . . 1 . . . . . . . .
#>  [8,] . . . 1 . 1 1 . 1 . . . . . . . . .
#>  [9,] . . 1 . 2 1 . 2 2 . . 1 1 . 1 . . 1
#> [10,] . 1 . . . . . . . . . 1 . . . . . .
#> [11,] . . . . 2 . . . . . 1 . . 2 . . . .
#> [12,] . . . 1 . . 1 . . . . . 1 . . . . 1
#> [13,] . . . . . . . . . . . 1 . 1 . 1 . .
#> [14,] . 1 . 1 . . . 1 . . . 1 . . . . . .
#> [15,] . . . . . . 1 . . . . . . . . . 1 .
#> [16,] . 1 1 . . . . . . . . 1 . . 1 . . .
#> [17,] . 1 . . . . . . . . . 1 . 1 . . 2 1
#> [18,] . . . 1 . . . 1 . 1 . . . . 1 1 . .
#> [19,] 1 1 . . . 1 1 . 1 1 . 1 . 1 . . 1 .
#> [20,] 1 . . . . . . . . . . . . . . . . 1
## same as cross-tabulation of the two index vectors:
stopifnot(all(c12 - unclass(table(IM1@perm, IM2@perm)) == 0))

# 3 observations, 4 implied values, first does not occur in sample:
as(2:4, "indMatrix")
#> 3 x 4 sparse Matrix of class "indMatrix"
#>             
#> [1,] . | . .
#> [2,] . . | .
#> [3,] . . . |
# 3 observations, 5 values, first and last do not occur in sample:
as(list(2:4, 5), "indMatrix")
#> 3 x 5 sparse Matrix of class "indMatrix"
#>               
#> [1,] . | . . .
#> [2,] . . | . .
#> [3,] . . . | .

as(sm1, "nMatrix")
#> 9 x 3 sparse Matrix of class "ngRMatrix"
#>            
#>  [1,] . | .
#>  [2,] . | .
#>  [3,] . | .
#>  [4,] . . |
#>  [5,] . . |
#>  [6,] . . |
#>  [7,] | . .
#>  [8,] | . .
#>  [9,] | . .
s10[1:7, 1:4] # gives an "ngTMatrix" (most economic!)
#> 7 x 4 sparse Matrix of class "ngCMatrix"
#>             
#> [1,] . | . .
#> [2,] . . . .
#> [3,] . . . .
#> [4,] . . . .
#> [5,] . . . .
#> [6,] . | . .
#> [7,] . . . .
s10[1:4, ]  # preserves "indMatrix"-class
#> 4 x 10 sparse Matrix of class "indMatrix"
#>                         
#> [1,] . | . . . . . . . .
#> [2,] . . . . . . . . . |
#> [3,] . . . . . . | . . .
#> [4,] . . . . | . . . . .

I1 <- as(c(5:1,6:4,7:3), "indMatrix")
I2 <- as(7:1, "pMatrix")
(I12 <- rbind(I1, I2))
#> 20 x 7 sparse Matrix of class "indMatrix"
#>                    
#>  [1,] . . . . | . .
#>  [2,] . . . | . . .
#>  [3,] . . | . . . .
#>  [4,] . | . . . . .
#>  [5,] | . . . . . .
#>  [6,] . . . . . | .
#>  [7,] . . . . | . .
#>  [8,] . . . | . . .
#>  [9,] . . . . . . |
#> [10,] . . . . . | .
#> [11,] . . . . | . .
#> [12,] . . . | . . .
#> [13,] . . | . . . .
#> [14,] . . . . . . |
#> [15,] . . . . . | .
#> [16,] . . . . | . .
#> [17,] . . . | . . .
#> [18,] . . | . . . .
#> [19,] . | . . . . .
#> [20,] | . . . . . .
stopifnot(is(I12, "indMatrix"),
          identical(I12, rbind(I1, I2)),
    colSums(I12) == c(2L,2:4,4:2))