[BioC] Gating with an ellipsoidGate

Wed Sep 15 13:24:06 CEST 2010

Dear Josef,

Many thanks, the speadsheet was most helpful!

I did try creating a covariance matrix from separate scaling and 
rotation matrices (below), but this didn't work, I suppose because they 
are too simple--operating in Euclidean space rather than Mahalanobis space.

cov.matrix.simple <- function (a, b, angle) {
   theta <- angle * (pi/180)

   R <- matrix(c(cos(theta), -sin(theta), sin(theta), cos(theta)), 
byrow=TRUE, ncol=2)
   S <- matrix(c(a, 0, 0, b), byrow=TRUE, ncol=2)

   S%*%R
}

Using your spreadsheet, I've created a macro to create a covariance 
matrix in R from a, b and an angle, which I've verified to match the 
results shown in your spreadsheet:

cov.matrix <- function (a, b, angle) {
   theta <- angle * (pi/180)

   c1 <- ((cos(theta)^2)/a^2) + ((sin(theta)^2)/b^2)
   c2 <- sin(theta) * cos(theta) * ((1/a^2) - (1/b^2))
   c3 <- ((sin(theta)^2)/a^2) + ((cos(theta)^2)/b^2)

   m1 <- matrix(c(c1, c2, c2, c3), byrow=TRUE, ncol=2)
   m2 <- solve(m1)

   m2
}

This appears to work well when rotation is not involved.  However, as 
soon as I try to add rotation, it blows up:

library(flowCore)
library(flowViz)
library(flowUtils)

cov.matrix <- function (a, b, angle) {
   theta <- angle * (pi/180)

   c1 <- ((cos(theta)^2)/a^2) + ((sin(theta)^2)/b^2)
   c2 <- sin(theta) * cos(theta) * ((1/a^2) - (1/b^2))
   c3 <- ((sin(theta)^2)/a^2) + ((cos(theta)^2)/b^2)

   m1 <- matrix(c(c1, c2, c2, c3), byrow=TRUE, ncol=2)
   m2 <- solve(m1)

   m2
}

d <- read.FCS("0isotype.fcs", alter.names=TRUE)
d <- transform(d, `SS.Log` = log10(`SS.Log`))

# No rotation
cov <- cov.matrix(20000, 0.4, 0)
colnames(cov) <- c("FS.Lin", "SS.Log")
rownames(cov) <- c("FS.Lin", "SS.Log")
mean <- c("FS.Lin"=40000, "SS.Log"=2.8)
cells <- ellipsoidGate(filterId="CellGate", .gate=cov, mean=mean)
print(cov)
pdf("test.pdf", width=8, height=8, pointsize=12)
print(xyplot(`SS.Log` ~ `FS.Lin`, d, filter=cells, xlab="FS", 
ylab=expression(log[10]~(SS))))
dev.off()

# Now repeat with a 20 degree rotation...
cov <- cov.matrix(20000, 0.4, 20)
colnames(cov) <- c("FS.Lin", "SS.Log")
rownames(cov) <- c("FS.Lin", "SS.Log")
mean <- c("FS.Lin"=40000, "SS.Log"=2.8)
cells <- ellipsoidGate(filterId="CellGate", .gate=cov, mean=mean)
print(cov)

pdf("test-rotate.pdf", width=8, height=8, pointsize=12)
print(xyplot(`SS.Log` ~ `FS.Lin`, d, filter=cells, xlab="FS", 
ylab=expression(log[10]~(SS))))
dev.off()

I've put the source data, script and results at:
http://www-users.york.ac.uk/~rl522/flowcore-test/

http://www-users.york.ac.uk/~rl522/flowcore-test/test.pdf
http://www-users.york.ac.uk/~rl522/flowcore-test/test-rotate.pdf

You can see that while an ellipse is drawn correctly in the first 
instance, when we add a 20 degree rotation, it's completely screwed. 
I'm not sure if this is due to the extreme differences in the x and y 
dimensions, which alter the length of a and b, or for some other reason.

Many thanks,
Roger

On 14/09/2010 18:00, Josef Spidlen wrote:
> Hi Roger,
> I am not the author of flowCore but I believe the ellipsoid gate
> constructor takes covariance matrices since this is used in the
> Gating-ML specification. Mathematically, this is a nice way to specify
> ellipsoids in multidimensional space.
>
> I have posted a simple spreadsheet that allows you to convert the
> representation of an ellipse by half-axes, rotation and centre point to
> covariance matrix, mean (= centre point), and distance square (=1 for
> this case). Please follow this link for download:
> http://sourceforge.net/projects/flowcyt/files/Gating-ML/Gating-ML%201.5/EllipseCalculations.xls/download
>
> I hope this addresses your issues; let me know if further help is needed.
>
> Unfortunately, I'll have to leave it up to the flowCore/flowViz authors
> to get back to you with your other questions.
>
> Cheers,
> Josef
>
>> Message: 6
>> Date: Mon, 13 Sep 2010 15:51:55 +0100
>> From: Roger Leigh<rleigh at codelibre.net>
>> To: bioconductor at stat.math.ethz.ch
>> Subject: [BioC] FlowCore/FlowViz issues
>> Message-ID:<20100913145155.GW6128 at codelibre.net>
>> Content-Type: text/plain; charset=utf-8
>>
>> Hi,
>>
>> [Apologies if this is delivered twice; my initial mail didn't
>> appear to be accepted for several hours, and I didn't see any
>> failure message; does the list object to GPG signatures?]
>>
>> I've just started to use FlowCore/FlowViz to analyse some of my
>> flow cytometry data, and ran into a few problems. I'm hoping
>> that you might be able to point me in the right direction!
>>
>> I've been very pleased with it so far, and have got some nice
>> plots and stats out of it, but I'm sure I'm doing some things
>> very inefficiently and/or incorrectly!
>>
>> [I'm using R version 2.11.1 (2010-05-31) on x86_64-pc-linux-gnu
>> (Debian GNU/Linux) with current Bioconductor packages)]
>>
>>
>> 1) Gating with an ellipsoidGate
>>
>> cov<- matrix(c(400000000, 0, 0, 0.08), ncol=2,
>> dimnames=list(c("FS.Lin", "SS.Log"), c("FS.Lin", "SS.Log")))
>> mean<- c("FS.Lin"=32000, "SS.Log"=2.8)
>> cells<- ellipsoidGate(filterId="CellGate", .gate=cov, mean=mean)
>>
>> I want to select my cells using an ellipsoid gate on forward- and side-
>> scatter plots. In the above situation, they lie in a region where
>> FS=32000?10000 and log??(SS)=2.8?0.5. However, the values in the
>> covariance matrix don't match the dimensions; is there any explanation
>> regarding how to construct a covariance matrix from the actual
>> dimension I want (I got the above by trial and error until it fitted
>> nicely--I'm afraid I know little about these matrices).
>>
>> In some plots I'd also like to rotate the ellipse, but I'm not sure
>> how to put this into the matrix, if that's the way to do things.
>> Is this possible?
>>
>> Is there an alternative constructor to create a gate from real
>> dimensions?
>>
>>
>>
>> Many thanks for all your help,
>> Roger
>>
>> --
>> .''`. Roger Leigh
>> : :' : Debian GNU/Linux http://people.debian.org/~rleigh/
>> `. `' Printing on GNU/Linux? http://gutenprint.sourceforge.net/
>> `- GPG Public Key: 0x25BFB848 Please GPG sign your mail.
>>
>>
>