# [Rd] matrix bands

From: Jeremy David Silver <jds_at_dmu.dk>
Date: Fri, 26 Aug 2011 13:23:43 +0200

Regards,
Jeremy Silver

```###############################################

```

## less clear formulation - more efficient band <- function(x,n = 0){

dx <- dim(x)
if(length(dx) != 2L)
stop("only matrix bands can be accessed")

max.dx <- max(dx)
n <- as.integer(n)

ij <- cbind(i = seq(1,max.dx) - n,

j = seq(1,max.dx))
ij <- ij[1 <= ij[,1] & ij[,1] <= dx[1] & 1 <= ij[,2] & ij[,2] <= dx[2],,drop=FALSE]

if(nrow(ij) == 0)
stop('cannot access this matrix band')

x[ij]
}

dx <- dim(x)
if(length(dx) != 2L)
stop("only matrix bands can be replaced")

max.dx <- max(dx)
n <- as.integer(n)

ij <- cbind(i = seq(1,max.dx) - n,

j = seq(1,max.dx))
ij <- ij[1 <= ij[,1] & ij[,1] <= dx[1] & 1 <= ij[,2] & ij[,2] <= dx[2],,drop=FALSE]

if(nrow(ij) == 0)
stop('cannot replace this matrix band')

x[ij] <- value

x
}

x[col(x) - row(x) == as.integer(n)]
}

## define a test matrix
> A <- matrix(rnorm(12),3,4)
> A

[,1] [,2] [,3] [,4]
[1,] -1.5560200 0.6452762 1.072565 0.1923451
[2,] 0.7940685 1.2441817 1.699486 -0.2998814

[3,] -0.7762252 -0.4824173 -0.981055 -0.9265627

## access some of the bands

> band(A,1)
[1] 0.6452762 1.6994858 -0.9265627

> band(A,-2)
[1] -0.7762252

> band(A,2)
[1] 1.0725649 -0.2998814

## set one of the bands

> band(A,2) <- 2:1
> A

[,1] [,2] [,3] [,4]
[1,] -1.5560200 0.6452762 2.000000 0.1923451
[2,] 0.7940685 1.2441817 1.699486 1.0000000

[3,] -0.7762252 -0.4824173 -0.981055 -0.9265627

> A <- matrix(1:10)
> A

[,1]

```  [1,]    1
[2,]    2
[3,]    3
[4,]    4
[5,]    5
[6,]    6
[7,]    7
[8,]    8
[9,]    9
```

[10,] 10

> band(A,0)
[1] 1

> band(A,1)
Error in band(A, 1) : cannot access this matrix band  > band(A,-1)
[1] 2

> band(A,-5)
[1] 6

> for(i in -2:3){print(band(A,i));print(band2(A,i))}
[1] -0.7762252
[1] -0.7762252
[1] 0.7940685 -0.4824173
[1] 0.7940685 -0.4824173
[1] -1.556020 1.244182 -0.981055
[1] -1.556020 1.244182 -0.981055

[1] 0.6452762 1.6994858 -0.9265627
[1] 0.6452762 1.6994858 -0.9265627
[1] 2 1
[1] 2 1
[1] 0.1923451
[1] 0.1923451

``` > N <- 1e4
> M <- matrix(0,N,N)
> system.time(band(M,2))
```

user system elapsed
0.005 0.003 0.007
> system.time(band2(M,2))

user system elapsed
18.509 2.121 20.754

R-devel_at_r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-devel Received on Fri 26 Aug 2011 - 11:41:24 GMT

This quarter's messages: by month, or sorted: [ by date ] [ by thread ] [ by subject ] [ by author ]

Archive maintained by Robert King, hosted by the discipline of statistics at the University of Newcastle, Australia.
Archive generated by hypermail 2.2.0, at Fri 26 Aug 2011 - 13:10:23 GMT.

Mailing list information is available at https://stat.ethz.ch/mailman/listinfo/r-devel. Please read the posting guide before posting to the list.