[R] iplots Java

[Michael Kubovy]

I run R 2.5.1 on a Mac, and use JGR as the front-end. When I  
performed an update.packages(), I believe it unpdated a component of  
JGR. I then quit and tried to relaunch JGR. It wouldn't launch.  
Instead it opened a panel that says: "Cannot find iplot Java classes.  
Please make sure that the latest iplots R package is correctly  
installed." I would appreeciate hearing of strategies for solving the  
problem.

_
Professor Michael Kubovy
University of Virginia
Department of Psychology
USPS: P.O.Box 400400Charlottesville, VA 22904-4400
Parcels:Room 102Gilmer Hall
 McCormick RoadCharlottesville, VA 22903
Office:B011+1-434-982-4729
Lab:B019+1-434-982-4751
Fax:+1-434-982-4766
WWW:http://www.people.virginia.edu/~mk9y/

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

Re: [R] Artifacts in pdf() of image() (w/o comments)

[Michael Kubovy]

Hi Duncan,

I was trying to learn to remove the artifacts by setting the  
parameters of the image so that anti-aliasing wouldn't produce them  
(analogous to making sure that the screens of two halftone screen  
process images are in register before combining them so as to avoid  
Moiré patterns). So, I wasn't complaining; I was looking for  
instruction. Since the image was computed with R, it's reasonable for  
me to ask R experts to help me out of what I thought was caused by my  
ineptitude. (Regarding your question, "If it wasn't a bug, why did it  
bother you so much?"---there are many things that bother me that are  
not bugs.)

Moreover, I'm not opposed to complaining to Apple, once I have been  
assured that I'm not reporting a bug where there's none.

On Aug 13, 2007, at 9:56 AM, Duncan Murdoch wrote:

> On 8/13/2007 11:43 AM, Michael Kubovy wrote:
>> But is it a bug? Can a program anti-alias text and line drawings  
>> and  not bitmaps?
>
> Anti-aliasing is the removal of artifacts caused by displaying an  
> image on a low-resolution bitmapped display.  Introducing artifacts  
> is a bug.
>
> If it wasn't a bug, why did it bother you so much?  And why do you  
> think it's reasonable to complain about it on R-help, but not to  
> complain about it to Apple, who are clearly responsible for it?
>
> Duncan Murdoch
>
>> On Aug 13, 2007, at 9:30 AM, Duncan Murdoch wrote:
>>> On 8/13/2007 11:07 AM, Michael Kubovy wrote:
>>>> Dear Friends,
>>>> Thanks for your input.
>>>> FYI: Preview doesn't show PDF aliasing in the image I produced  
>>>> if  I  uncheck the "Anti-alias text and line art" box under the  
>>>> PDF  tab in  Preferences. So I'm not yet ready to drop Preview  
>>>> from my  toolbox.
>>>
>>> An alternative to dropping Preview is to report the bug in it to   
>>> Apple.  Apple has an online bug reporting web page somewhere; I   
>>> haven't found them as helpful as R-help, but your mileage may vary.
>>>
>>> Duncan Murdoch
>> _
>> Professor Michael Kubovy
>> University of Virginia
>> Department of Psychology
>> USPS: P.O.Box 400400Charlottesville, VA 22904-4400
>> Parcels:Room 102Gilmer Hall
>>  McCormick RoadCharlottesville, VA 22903
>> Office:B011+1-434-982-4729
>> Lab:B019+1-434-982-4751
>> Fax:+1-434-982-4766
>> WWW:http://www.people.virginia.edu/~mk9y/

_
Professor Michael Kubovy
University of Virginia
Department of Psychology
USPS: P.O.Box 400400Charlottesville, VA 22904-4400
Parcels:Room 102Gilmer Hall
 McCormick RoadCharlottesville, VA 22903
Office:B011+1-434-982-4729
Lab:B019+1-434-982-4751
Fax:+1-434-982-4766
WWW:http://www.people.virginia.edu/~mk9y/

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

Re: [R] Artifacts in pdf() of image() (w/o comments)

[Michael Kubovy]

But is it a bug? Can a program anti-alias text and line drawings and  
not bitmaps?

On Aug 13, 2007, at 9:30 AM, Duncan Murdoch wrote:

> On 8/13/2007 11:07 AM, Michael Kubovy wrote:
>> Dear Friends,
>> Thanks for your input.
>> FYI: Preview doesn't show PDF aliasing in the image I produced if  
>> I  uncheck the "Anti-alias text and line art" box under the PDF  
>> tab in  Preferences. So I'm not yet ready to drop Preview from my  
>> toolbox.
>
> An alternative to dropping Preview is to report the bug in it to  
> Apple.  Apple has an online bug reporting web page somewhere; I  
> haven't found them as helpful as R-help, but your mileage may vary.
>
> Duncan Murdoch

_
Professor Michael Kubovy
University of Virginia
Department of Psychology
USPS: P.O.Box 400400Charlottesville, VA 22904-4400
Parcels:Room 102Gilmer Hall
 McCormick RoadCharlottesville, VA 22903
Office:B011+1-434-982-4729
Lab:B019+1-434-982-4751
Fax:+1-434-982-4766
WWW:http://www.people.virginia.edu/~mk9y/

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

Re: [R] Artifacts in pdf() of image() (w/o comments)

[Michael Kubovy]

Dear Friends,

Thanks for your input.

FYI: Preview doesn't show PDF aliasing in the image I produced if I  
uncheck the "Anti-alias text and line art" box under the PDF tab in  
Preferences. So I'm not yet ready to drop Preview from my toolbox.

MK

On Aug 13, 2007, at 12:02 AM, Mark Wardle wrote:

> It may be worth outputting postscript and converting to PDF from
> there. Although Preview can do this, it may be worth looking at
> Ghostscript which *may* not have simila problems.
>
> I have also had PDFs which have displayed well in Preview and open
> source tools and have been garbled in Adobe Acrobat, so the problems
> aren't limited to Preview.
>
> Best wishes,
>
> Mark
>
> On 13/08/07, Duncan Murdoch <[EMAIL PROTECTED]> wrote:
>
>> Michael Kubovy wrote:
>>> On Aug 12, 2007, at 6:24 AM, Duncan Murdoch wrote:
>>>
>>>> Michael Kubovy wrote:
>>>>
>>>>> Dear r-helpers,
>>>>>
>>>>> In my previous message there were comments in the code that may
>>>>> have  made cutting and pasting awkward. Here it is w/o them.
>>>>>
>>>>> I have two questions:
>>>>>
>>>>> (1) The following produces a pdf with artifacts. How do I prevent
>>>>> them?
>>>>>
>>>>>
>>>> What artifacts do you see?  It looks like a smoothly varying field
>>>> when produced by R 2.5.1 and viewed in Acrobat Reader 6.0 on  
>>>> Windows.
>>>>
>>>> Duncan Murdoch
>>>>
>>>>> require(grDevices)
>>>>> imSize <- 200
>>>>> lambda <- 10
>>>>> theta <- 15
>>>>> sigma <- 40
>>>>> x <- 1:imSize
>>>>> x0 <- x / imSize -.5
>>>>> freq = imSize/lambda
>>>>> xf = x0 * freq * 2 * pi
>>>>> f <- function(x, y){r <- -((x^2 + y^2)/(sigma ^2)); exp(r)}
>>>>> z <- outer(xf, xf, f)
>>>>> f1 <- function(x, y){cos(.1 * x)}
>>>>> z1 <- outer(xf, xf, f1)
>>>>> pdf('gabor.pdf')
>>>>> image(xf, xf, z * z1, col = gray(250:1000/1000),
>>>>> xlab = '', ylab =  '', bty = 'n', axes = FALSE, asp = 1)
>>>>> dev.off()
>>>>>
>>>
>>> I'm working on a Mac. You're right, Acrobat 6.05 renders the figure
>>> nicely, but when it's included in a LaTeX-produced pdf or viewed  
>>> with
>>> the Mac Preview program, a grid of fine white lines is superimposed
>>> on the figure. So I believe that it's a matter of aliasing, which I
>>> might be able to prevent by adjusting the parameters of the figures.
>>> I just don't know enough to figure this out, and would appreciate
>>> guidance.
>> I see the artifacts in Preview on a Mac too.  So it looks to me  
>> like a
>> Mac bug.
>>
>> Preview is actually pretty poor at graphics display; see
>> <http://www.geuz.org/pipermail/gl2ps/2007/000223.html>.
>>
>> My only suggestion is not to use Preview.
>>
>> Duncan Murdoch
> -- 
> Dr. Mark Wardle
> Clinical research fellow and specialist registrar, Neurology
> Cardiff, UK

_
Professor Michael Kubovy
University of Virginia
Department of Psychology
USPS: P.O.Box 400400Charlottesville, VA 22904-4400
Parcels:Room 102Gilmer Hall
 McCormick RoadCharlottesville, VA 22903
Office:B011+1-434-982-4729
Lab:B019+1-434-982-4751
Fax:+1-434-982-4766
WWW:http://www.people.virginia.edu/~mk9y/

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

Re: [R] Artifacts in pdf() of image() (w/o comments)

[Michael Kubovy]

On Aug 12, 2007, at 6:24 AM, Duncan Murdoch wrote:

> Michael Kubovy wrote:
>> Dear r-helpers,
>>
>> In my previous message there were comments in the code that may  
>> have  made cutting and pasting awkward. Here it is w/o them.
>>
>> I have two questions:
>>
>> (1) The following produces a pdf with artifacts. How do I prevent  
>> them?
>>
>
> What artifacts do you see?  It looks like a smoothly varying field  
> when produced by R 2.5.1 and viewed in Acrobat Reader 6.0 on Windows.
>
> Duncan Murdoch
>> require(grDevices)
>> imSize <- 200
>> lambda <- 10
>> theta <- 15
>> sigma <- 40
>> x <- 1:imSize
>> x0 <- x / imSize -.5
>> freq = imSize/lambda
>> xf = x0 * freq * 2 * pi
>> f <- function(x, y){r <- -((x^2 + y^2)/(sigma ^2)); exp(r)}
>> z <- outer(xf, xf, f)
>> f1 <- function(x, y){cos(.1 * x)}
>> z1 <- outer(xf, xf, f1)
>> pdf('gabor.pdf')
>> image(xf, xf, z * z1, col = gray(250:1000/1000),
>> xlab = '', ylab =  '', bty = 'n', axes = FALSE, asp = 1)
>> dev.off()

I'm working on a Mac. You're right, Acrobat 6.05 renders the figure  
nicely, but when it's included in a LaTeX-produced pdf or viewed with  
the Mac Preview program, a grid of fine white lines is superimposed  
on the figure. So I believe that it's a matter of aliasing, which I  
might be able to prevent by adjusting the parameters of the figures.  
I just don't know enough to figure this out, and would appreciate  
guidance.

_
Professor Michael Kubovy
University of Virginia
Department of Psychology
USPS: P.O.Box 400400Charlottesville, VA 22904-4400
Parcels:Room 102Gilmer Hall
 McCormick RoadCharlottesville, VA 22903
Office:B011+1-434-982-4729
Lab:B019+1-434-982-4751
Fax:+1-434-982-4766
WWW:http://www.people.virginia.edu/~mk9y/

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

[R] Artifacts in pdf() of image() (w/o comments)

[Michael Kubovy]

Dear r-helpers,

In my previous message there were comments in the code that may have  
made cutting and pasting awkward. Here it is w/o them.

I have two questions:

(1) The following produces a pdf with artifacts. How do I prevent them?

require(grDevices)
imSize <- 200
lambda <- 10
theta <- 15
sigma <- 40
x <- 1:imSize
x0 <- x / imSize -.5
freq = imSize/lambda
xf = x0 * freq * 2 * pi
f <- function(x, y){r <- -((x^2 + y^2)/(sigma ^2)); exp(r)}
z <- outer(xf, xf, f)
f1 <- function(x, y){cos(.1 * x)}
z1 <- outer(xf, xf, f1)
pdf('gabor.pdf')
image(xf, xf, z * z1, col = gray(250:1000/1000),
xlab = '', ylab =  '', bty = 'n', axes = FALSE, asp = 1)
dev.off()

(2) I would like the output to be clipped to a circle, i.e., anything
outside the circle tangent to the sides of the square should be
transparent. How can I do that?

_
Professor Michael Kubovy
University of Virginia
Department of Psychology
USPS: P.O.Box 400400Charlottesville, VA 22904-4400
Parcels:Room 102Gilmer Hall
  McCormick RoadCharlottesville, VA 22903
Office:B011+1-434-982-4729
Lab:B019+1-434-982-4751
Fax:+1-434-982-4766
WWW:http://www.people.virginia.edu/~mk9y/

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting- 
guide.html
and provide commented, minimal, self-contained, reproducible code.

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

[R] Artifacts in pdf() of image()

[Michael Kubovy]

Dear r-helpers,

I have two questions:

(1) The following produces a pdf with artifacts. How do I prevent them?

require(grDevices)
imSize <- 200
lambda <- 10
theta <- 15
sigma <- 40
x <- 1:imSize
x0 <- x / imSize -.5
freq = imSize/lambda# compute frequency from  
wavelength
xf = x0 * freq * 2 * pi  # convert X to radians: 0 ->  
( 2*pi * frequency)
f <- function(x, y){r <- -((x^2 + y^2)/(sigma ^2)); exp(r)}
z <- outer(xf, xf, f)
f1 <- function(x, y){cos(.1 * x)}
z1 <- outer(xf, xf, f1)
pdf('gabor.pdf')
image(xf, xf, z * z1, col = gray(250:1000/1000), xlab = '', ylab =  
'', bty = 'n',
 axes = FALSE, asp = 1)
dev.off()

(2) I would like the output to be clipped to a circle, i.e., anything  
outside the circle tangent to the sides of the square should be  
transparent. How can I do that?

_
Professor Michael Kubovy
University of Virginia
Department of Psychology
USPS: P.O.Box 400400Charlottesville, VA 22904-4400
Parcels:Room 102Gilmer Hall
 McCormick RoadCharlottesville, VA 22903
Office:B011+1-434-982-4729
Lab:B019+1-434-982-4751
Fax:+1-434-982-4766
WWW:http://www.people.virginia.edu/~mk9y/

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

Re: [R] Standard errors of the predicted values from a lme (or lmer)-object

[Michael Kubovy]

On Jun 1, 2007, at 6:08 AM, Dieter Menne wrote:

> Fränzi Korner  oikostat.ch> writes:
>
>> how do I obtain standard errors of the predicted values from a lme  
>> (or
>> lmer)-object?
>
> Not totally clear what you mean. intervals(lmeresult) gives the  
> confidence
> intervals for the coefficients. Otherwise, you can do some  
> calculations with
> residuals(lmeresult). Most useful for diagnostic purposes is plot 
> (lmeresult).

Perhaps
?estimable
in the gmodels package
_____
Professor Michael Kubovy
University of Virginia
Department of Psychology
USPS: P.O.Box 400400Charlottesville, VA 22904-4400
Parcels:Room 102Gilmer Hall
 McCormick RoadCharlottesville, VA 22903
Office:B011+1-434-982-4729
Lab:B019+1-434-982-4751
Fax:+1-434-982-4766
WWW:http://www.people.virginia.edu/~mk9y/

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

Re: [R] Selecting complementary colours

[Michael Kubovy]

Hi All,

Complementary colors will not insure legibility. If the text and the  
background are equiluminant, visibility will be severely impaired.

On May 21, 2007, at 8:22 PM, John Fox wrote:

> Dear Achim,
>
> As I mentioned in my response to Deepayan's suggestion, I didn't  
> specify the
> original problem clearly: The object is to get contrasting colours,  
> so that
> when one is plotted over the other, the two will be readily  
> distinguishable.
> Your suggestions don't do this for neutral colours:
>
>> x <- "#88"
>> y_hcl <- as(hex2RGB(x), "polarLUV")
>> [EMAIL PROTECTED], "H"] <- [EMAIL PROTECTED], "H"] + 180
>> hex(y_hcl)
> [1] "#88"
>
>> y_hsv <- as(hex2RGB(x), "HSV")
>> [EMAIL PROTECTED], "H"] <- [EMAIL PROTECTED], "H"] + 180
>> hex(y_hsv)
> [1] "#88"
>
> Thank you for trying to help.
>
> John
>
> 
> John Fox, Professor
> Department of Sociology
> McMaster University
> Hamilton, Ontario
> Canada L8S 4M4
> 905-525-9140x23604
> http://socserv.mcmaster.ca/jfox
> 
>
>> -Original Message-
>> From: Achim Zeileis [mailto:[EMAIL PROTECTED]
>> Sent: Monday, May 21, 2007 7:07 PM
>> To: John Fox
>> Cc: r-help@stat.math.ethz.ch
>> Subject: Re: [R] Selecting complementary colours
>>
>> On Mon, 21 May 2007, John Fox wrote:
>>
>>> Dear r-helpers,
>>>
>>> I wonder whether, given the "#rrggbb" representation of a colour,
>>> there is a simple way to select the complementary colour,
>> also expressed as a "#rrggbb"
>>> string.
>>
>> Is the complementary color uniquely defined? My understanding
>> is that you can take opposite colors on a color wheel, but
>> there are of course various color wheels available. With
>> "colorspace" you can experiment with this,
>> e.g.:
>>   x <- "#81A9D0"
>>   y_hcl <- as(hex2RGB(x), "polarLUV")
>>   [EMAIL PROTECTED], "H"] <- [EMAIL PROTECTED], "H"] + 180
>>   y_hcl <- hex(y_hcl)
>> which is a bit more balanced than
>>   y_hsv <- as(hex2RGB(x), "HSV")
>>   [EMAIL PROTECTED], "H"] <- [EMAIL PROTECTED], "H"] + 180
>>   y_hsv <- hex(y_hsv)
>>
>> hth,
>> Z
>>
>>
>>
>
> __
> R-help@stat.math.ethz.ch mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide http://www.R-project.org/posting- 
> guide.html
> and provide commented, minimal, self-contained, reproducible code.

_
Professor Michael Kubovy
University of Virginia
Department of Psychology
USPS: P.O.Box 400400Charlottesville, VA 22904-4400
Parcels:Room 102Gilmer Hall
 McCormick RoadCharlottesville, VA 22903
Office:B011+1-434-982-4729
Lab:B019+1-434-982-4751
Fax:+1-434-982-4766
WWW:http://www.people.virginia.edu/~mk9y/

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

Re: [R] title line

[Michael Kubovy]

On May 21, 2007, at 8:21 AM, elyakhlifi mustapha wrote:

> I know how to do to write down the columns in a data frame but I  
> also wanna write down the lines .
> Do you know how could I do this?

I'm not sure what you mean by "write down the lines". If you mean  
that you want to know what the row names are, then an answer is:

df <- data.frame(x1 = rnorm(3), x2 = rnorm(3))
names(df) # just to show that I gave the columns names while defining df
rownames(df)

If you mean that you want to assign names to the rows, then

df <- data.frame(rnorm(3), rnorm(3))
names(df) <- c('y1', 'y2') # not required
rownames(df) <- c('a', 'b', 'c')

In one commandt:

df <- data.frame(x1 = rnorm(3), x2 = rnorm(3), row.names = c('a',  
'b', 'c'))

More detailed:

df <- data.frame(x1 = rnorm(3), x2 = rnorm(3), row.names = c('a',  
'b', 'c'), check.rows = TRUE, check.names = TRUE)


More importantly you could have found all this out by typing ?data.frame


I have seen several posts from you that seem to be rather elementary  
questions, sometimes not using the correct terminology ('write down'  
and 'lines of a data frame' are not standard). Have you taken the  
time to go through any of the many introductions to R available on line?

Also, you're not asking your questions well (this comment is not  
about knowing English well):
You would have been better off asking:
"I have a data frame
df <- data.frame(x1 = rnorm(3), x2 = rnorm(3))
How do I ? "
_
Professor Michael Kubovy
University of Virginia
Department of Psychology
USPS: P.O.Box 400400Charlottesville, VA 22904-4400
Parcels:Room 102Gilmer Hall
 McCormick RoadCharlottesville, VA 22903
Office:B011+1-434-982-4729
Lab:B019+1-434-982-4751
Fax:+1-434-982-4766
WWW:http://www.people.virginia.edu/~mk9y/

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

[R] Separate x-axis ticks for panels in xyplot()?

[Michael Kubovy]

Dear r-helpers,

In an xyplot I have

xyplot(y ~ x | c, groups = g,
scales = list(
x = list(
at = v1,
labels = c('A', 'B', 'C', 'D'),
rot = 45
)
)
)

g consists of two groups. How do I define a different set of x ticks  
and labels for the second panel? I'm looking for something like

xyplot(y ~ x | c, groups = g,
scales = list(
x = list(ifelse(g == g1, {
at = v1,
labels = c('A', 'B', 'C', 'D'),
rot = 45},
{at = v2,
labels = c('E', 'F' 'G"),
    rot = 45}
)
)
)


_
Professor Michael Kubovy
University of Virginia
Department of Psychology
USPS: P.O.Box 400400Charlottesville, VA 22904-4400
Parcels:Room 102Gilmer Hall
 McCormick RoadCharlottesville, VA 22903
Office:B011+1-434-982-4729
Lab:B019+1-434-982-4751
Fax:+1-434-982-4766
WWW:http://www.people.virginia.edu/~mk9y/

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

Re: [R] xyplot() and controlling panel.polygon()

[Michael Kubovy]

I can't thank you enough! (I am also grateful to Sundar Dorai-Raj and  
Bert Gunter for their initial responses.)

On Apr 27, 2007, at 4:14 PM, Deepayan Sarkar wrote:

> On 4/27/07, Michael Kubovy <[EMAIL PROTECTED]> wrote:
>> Hi Deepayan,
>>
>> Your solution works, anf the polygon are drawn where I wanted them to
>> go. I thought that I could figure out how to gain control over the
>> colors of the four ensuing polygons (I'm trying to get two lighter
>> shades of the lines).
>>
>> I've tried, for example, to see if I could control the color of the
>> polyon outline, by adding border = 'red' to panel.polygon. That
>> didn't work. Does it work only in lpolygon()?
>
> No. This is a bug, and the border color can currently only be black or
> transparent (but it's easy to find and fix the bug; lpolygon is very
> short, so just write a replacement).
>
> However, panel.superpose recognizes and splits certain graphical
> parameters; border is not among them, but fill is. So you could do:
>
> my.panel.bands <-
>function(x, y, upper, lower,
> fill, col,
> subscripts, ..., font, fontface)
> {
>upper <- upper[subscripts]
>lower <- lower[subscripts]
>panel.polygon(c(x, rev(x)), c(upper, rev(lower)),
>  col = fill, border = FALSE,
>  ...)
> }
>
> xyplot(est ~ x | cond, group = grp, data = data, type = 'b',
>   col = c("#0080ff", "#ff00ff"),
>   fill = c("#bbddff", "#ffbbff"),
>   upper = data$upper,
>   lower = data$lower,
>   panel = function(x, y, ...){
>   panel.superpose(x, y, panel.groups = 'my.panel.bands', ...)
>   panel.xyplot(x, y, ...)
>   })
>
>
>> I often can figure things out on my own, but obviously there's
>> something fundamental that I'm not getting about inheritance and
>> passing in these sorts of objects. I've been trying to get it from
>> the help pages and from Murrell's book, but neither offers enough of
>> a cookbook for me to figure these things out. Is there something I
>> should have read?
>
> The concepts are all there in the help page (but it's often difficult
> to put them together). The main points are:
>
> (1) unrecognized arguments get passed on to the panel function as is
> (2) subscripts give indices of x, y, etc in the original data frame
>
> (the implication being that if you have another column from the
> original data frame, such as upper and lower in your example, indexing
> by subscripts will give you the matching subset). Other than that, the
> panel functions have to do their own work (and what they do should
> ideally be documented); nothing is enforced, so nothing is guaranteed.
>
> This sort of thing doesn't get used often enough for examples to be
> easily found. The following demo in lattice might be helpful.
>
> file.show(system.file("demo/intervals.R", package = "lattice"))
>
> -Deepayan

_
Professor Michael Kubovy
University of Virginia
Department of Psychology
USPS: P.O.Box 400400Charlottesville, VA 22904-4400
Parcels:Room 102Gilmer Hall
 McCormick RoadCharlottesville, VA 22903
Office:B011+1-434-982-4729
Lab:B019+1-434-982-4751
Fax:+1-434-982-4766
WWW:http://www.people.virginia.edu/~mk9y/

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

Re: [R] xyplot() and controlling panel.polygon()

[Michael Kubovy]

Hi Deepayan,

Your solution works, anf the polygon are drawn where I wanted them to  
go. I thought that I could figure out how to gain control over the  
colors of the four ensuing polygons (I'm trying to get two lighter  
shades of the lines).

I've tried, for example, to see if I could control the color of the  
polyon outline, by adding border = 'red' to panel.polygon. That  
didn't work. Does it work only in lpolygon()?

I often can figure things out on my own, but obviously there's  
something fundamental that I'm not getting about inheritance and  
passing in these sorts of objects. I've been trying to get it from  
the help pages and from Murrell's book, but neither offers enough of  
a cookbook for me to figure these things out. Is there something I  
should have read?

Thanks,
MIchael

On Apr 27, 2007, at 2:25 PM, Deepayan Sarkar wrote:

> On 4/27/07, Michael Kubovy <[EMAIL PROTECTED]> wrote:
>> Hi Deepayan,
>>
>> Thanks for your advice. This is moving along, however:
>> The following is drawing the same polygons in each panel. I'm trying
>> to get a different polygon (confidence band) for each group in each
>> panel. That's why I thought I would need to pass groups and
>> subscripts to the panel.groups
>>
>> est <- c(1:4, 3:6, 7, 9, 11, 13, 12, 15, 18, 21)
>> cond <- rep(c('a','b'), each = 8)
>> grp <- rep(c('I', 'II'), each = 4, 2)
>> x <- rep(c(.5, .7, .9, 1.1), 4)
>> upper <- est + 1
>> lower = est - 1
>> data <- data.frame(est = est, x = x, cond = cond, grp = grp, upper =
>> upper, lower = lower)
>>
>> my.panel.polygon <- function(..., font, fontface)
>> {
>> panel.polygon(...)
>> }
>>
>> xyplot(est ~ x | cond, group = grp, data = data, type = 'b',
>>  panel = function(x, y, ...){
>>  panel.superpose(c(x, rev(x)), c(upper, rev(lower)),
>>  panel.groups = 'my.panel.polygon', default.units =
>> 'native', ...)
>>  panel.xyplot(x, y, ...)
>>  }
>> )
>>
>> It's pretty clear that panel.superpose is not getting its x and y
>> values after they are split by group and panel.
>
> You are not even trying to do that; you have
>
>  panel.superpose(c(x, rev(x)), c(upper, rev(lower)), <...>
>
> so your x=c(x, rev(x)) is not the same length as x (and  
> subscripts), and your
> y = c(upper, rev(lower)) is not the same length as anything. Also,
> your upper and lower are being taken from the global env, not data
> (they happen to be the same, but since your data has them, I assume
> you want to use them).
>
> Perhaps you are looking for something like this:
>
>
> panel.bands <-
>function(x, y, upper, lower,
> subscripts, ..., font, fontface)
> {
>upper <- upper[subscripts]
>lower <- lower[subscripts]
>panel.polygon(c(x, rev(x)), c(upper, rev(lower)), ...)
> }
>
> xyplot(est ~ x | cond, group = grp, data = data, type = 'b',
>   upper = data$upper,
>   lower = data$lower,
>   panel = function(x, y, ...){
>   panel.superpose(x, y, panel.groups = 'panel.bands', ...)
>   panel.xyplot(x, y, ...)
>   })
>
> -Deepayan

_
Professor Michael Kubovy
University of Virginia
Department of Psychology
USPS: P.O.Box 400400Charlottesville, VA 22904-4400
Parcels:Room 102Gilmer Hall
 McCormick RoadCharlottesville, VA 22903
Office:B011+1-434-982-4729
Lab:B019+1-434-982-4751
Fax:+1-434-982-4766
WWW:http://www.people.virginia.edu/~mk9y/

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

Re: [R] xyplot() and controlling panel.polygon()

[Michael Kubovy]

Hi Deepayan,

Thanks for your advice. This is moving along, however:
The following is drawing the same polygons in each panel. I'm trying  
to get a different polygon (confidence band) for each group in each  
panel. That's why I thought I would need to pass groups and  
subscripts to the panel.groups

est <- c(1:4, 3:6, 7, 9, 11, 13, 12, 15, 18, 21)
cond <- rep(c('a','b'), each = 8)
grp <- rep(c('I', 'II'), each = 4, 2)
x <- rep(c(.5, .7, .9, 1.1), 4)
upper <- est + 1
lower = est - 1
data <- data.frame(est = est, x = x, cond = cond, grp = grp, upper =
upper, lower = lower)

my.panel.polygon <- function(..., font, fontface)
{
panel.polygon(...)
}

xyplot(est ~ x | cond, group = grp, data = data, type = 'b',
 panel = function(x, y, ...){
 panel.superpose(c(x, rev(x)), c(upper, rev(lower)),
 panel.groups = 'my.panel.polygon', default.units =  
'native', ...)
 panel.xyplot(x, y, ...)
 }
)

It's pretty clear that panel.superpose is not getting its x and y  
values after they are split by group and panel.

_
Professor Michael Kubovy
University of Virginia
Department of Psychology
USPS: P.O.Box 400400Charlottesville, VA 22904-4400
Parcels:Room 102Gilmer Hall
 McCormick RoadCharlottesville, VA 22903
Office:B011+1-434-982-4729
Lab:B019+1-434-982-4751
Fax:+1-434-982-4766
WWW:http://www.people.virginia.edu/~mk9y/

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

[R] Example of mcmcsamp() failing with lmer() output

[Michael Kubovy]

Hi,

I would appreciate help with the following model

<<1>>=
gunload <- read.table(hh('datasets/gunload.dat'), header = T)
gunload$method <- factor(gunload$method, labels = c('new', 'old'))
gunload$physique <- factor(gunload$group, labels = c('slight',  
'average', 'heavy'))
gunload$team9 <- factor(rep(1:9, each = 2))
@

This is a split plot design with randomly selected teams for each  
physique, which was crossed with two methods.

Following the Oats example in chap. 1 of MEMSS, and adapting to lmer:

<<2>>=
require(lme4)
gl.lmer <- lmer(rounds ~ method * physique + (1 | physique/team),  
data = gunload)
mcmcsamp(gl.lmer)
@
lmer runs as expected, but mcmcsamp() gives:
Error: Leading minor of order 1 in downdated X'X is not positive  
definite
Error in t(.Call(mer_MCMCsamp, object, saveb, n, trans, verbose,  
deviance)) :
error in evaluating the argument 'x' in selecting a method for  
function 't'

I tried:

<<3>>=
require(nlme)
gl.lme <- lme(rounds ~ method * physique, random = ~1 | physique/ 
team, data = gunload)
@
which gave all the expected output except that two p-values were NaN.

Which aspect of the data is making the simulation fail, and is there  
anything to be done about it?

_
Professor Michael Kubovy
University of Virginia
Department of Psychology
USPS: P.O.Box 400400Charlottesville, VA 22904-4400
Parcels:Room 102Gilmer Hall
 McCormick RoadCharlottesville, VA 22903
Office:B011+1-434-982-4729
Lab:B019+1-434-982-4751
Fax:+1-434-982-4766
WWW:http://www.people.virginia.edu/~mk9y/

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

Re: [R] xyplot() and controlling panel.polygon()

[Michael Kubovy]

Hi Sundar,

On Apr 26, 2007, at 10:03 PM, Deepayan Sarkar wrote:

> On 4/26/07, Michael Kubovy <[EMAIL PROTECTED]> wrote:
>> Hi Sundar,
>>
>> Thanks for taking a look at this
>>
>> On Apr 26, 2007, at 1:05 PM, Sundar Dorai-Raj wrote:
>>
>> >
>> >
>> > Michael Kubovy said the following on 4/26/2007 7:20 AM:
>> >> Dear R-helpers,
>> >> How do I tell panel.polygon what greoup and panel it applies to
>> >> whithin
>> >> xyplot(y ~ x | c, groups = g
>> >>  panel = function(x, y, groups, ...){
>> >>  panel.polygon(x = xpol[c, g], y = ypol[c, g],
>> >> default.units  = 'native')
>> >>  panel.xYplot(x, y, groups, ...)
>> >>  llines(x = c(1, 6),  y = c(-24.28, 35.941667),  
>> lwd = 2, lty
>> >> =  3, col = 4)
>> >>  }
>> >> x[c, g] and y[c, g] describe the polygon I want plotted for group
>> >> g  in panel c.
>> >> _
>> >
>> > I believe you can just do:
>> >
>> > panel = function(x, y, subscripts, groups, ...) {
>> >   ## note addition of subscripts to arguments
>> >   panel.superpose(xpol, ypol, subscripts, groups,
>> >   panel.groups = "panel.polygon",
>> >   default.units = "native", ...)
>> >   panel.xYplot(x, y, subscripts, groups, ...)
>> >   llines(x = c(1, 6),  y = c(-24.28, 35.941667),
>> > lwd = 2, lty = 3, col = 4)
>> > }
>> >
>> > Also, it would be easier to provide tested answers if you gave a
>> > reproducible example.
>>
>> Here is a reproducible example in which I tried (unsuccesfully) to
>> apply your suggestion:
>>
>> est <- c(1:4, 3:6, 7, 9, 11, 13, 12, 15, 18, 21)
>> cond <- rep(c('a','b'), each = 8)
>> grp <- rep(c('I', 'II'), each = 4, 2)
>> x <- rep(c(.5, .7, .9, 1.1), 4)
>> upper <- est + 1
>> lower = est - 1
>> data <- data.frame(est = est, x = x, cond = cond, grp = grp, upper =
>> upper, lower = lower)
>>
>> xyplot(est ~ x | cond, group = grp, data = data, type = 'b',
>>  panel = function(x, y, subscripts, groups, ...){
>>  panel.superpose(c(x, rev(x)), c(upper, rev(lower)),
>> subscripts, groups,
>>  panel.groups = 'panel.polygon', default.units =
>> 'native', ...)
>>  panel.xyplot(x, y, subscripts, groups, ...)
>>  }
>> )
>> # returns: "Error in validGP(list(...)) : Must specify only one of
>> 'font' and 'fontface'", which is incomprehensible to me.
>
> Yeah, it's happening because of this:
>
>> grid::gpar(font = 1, fontface = NULL)
> Error in validGP(list(...)) : Must specify only one of 'font' and  
> 'fontface'
>
> panel.polygon needs a fix. Replace it by this:
>
> my.panel.polygon <- function(..., font, fontface)
> {
>panel.polygon(...)
> }
>
> xyplot(est ~ x | cond, group = grp, data = data, type = 'b',
>   panel = function(x, y, ...) {
>   panel.superpose(c(x, rev(x)), c(upper, rev(lower)),
>   panel.groups = my.panel.polygon,
>   ...)
>   panel.xyplot(x, y, ...)
>   })
>
> You are also better off not having arguments you are not using (your
> code would have failed because of incorrect positional matching).

I'm not sure I understand.

I'm trying to add a filled band to show CIs around the data (which is  
why I tried to use xYplot, where 'filled bands' have bug).
So with these data:

est <- c(1:4, 3:6, 7, 9, 11, 13, 12, 15, 18, 21)
cond <- rep(c('a','b'), each = 8)
grp <- rep(c('I', 'II'), each = 4, 2)
x <- rep(c(.5, .7, .9, 1.1), 4)
upper <- est + 1
lower = est - 1
data <- data.frame(est = est, x = x, cond = cond, grp = grp, upper =
upper, lower = lower)

  I'm trying to decorate this with bands:

xyplot(est ~ x | cond, group = grp, data = data, type = 'b')

So I'm not sure why your proposal is failing:

my.panel.polygon <- function(..., font, fontface)
{
panel.polygon(...)
}

xYplot(est ~ x | cond, group = grp, data = data, type = 'b',
 panel = function(x, y, subscripts, groups, ...){
 panel.superpose(c(x, rev(x)), c(upper, rev(lower)),  
subscripts, groups,
 panel.groups = 'my.panel.polygon', default.units =  
'native', ...)
 panel.xYplot(x, y,,...)
 }
)

_
Professor Michael Kubovy
University of Virginia
Department of Psychology
USPS: P.O.Box 400400Charlottesville, VA 22904-4400
Parcels:Room 102Gilmer Hall
 McCormick RoadCharlottesville, VA 22903
Office:B011+1-434-982-4729
Lab:B019+1-434-982-4751
Fax:+1-434-982-4766
WWW:http://www.people.virginia.edu/~mk9y/

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

Re: [R] xyplot() and controlling panel.polygon()

[Michael Kubovy]

Hi Sundar,

Thanks for taking a look at this

On Apr 26, 2007, at 1:05 PM, Sundar Dorai-Raj wrote:

>
>
> Michael Kubovy said the following on 4/26/2007 7:20 AM:
>> Dear R-helpers,
>> How do I tell panel.polygon what greoup and panel it applies to  
>> whithin
>> xyplot(y ~ x | c, groups = g
>>  panel = function(x, y, groups, ...){
>>  panel.polygon(x = xpol[c, g], y = ypol[c, g],  
>> default.units  = 'native')
>>  panel.xYplot(x, y, groups, ...)
>>  llines(x = c(1, 6),  y = c(-24.28, 35.941667), lwd = 2, lty 
>>  
>> =  3, col = 4)
>>  }
>> x[c, g] and y[c, g] describe the polygon I want plotted for group  
>> g  in panel c.
>> _
>
> I believe you can just do:
>
> panel = function(x, y, subscripts, groups, ...) {
>   ## note addition of subscripts to arguments
>   panel.superpose(xpol, ypol, subscripts, groups,
>   panel.groups = "panel.polygon",
>   default.units = "native", ...)
>   panel.xYplot(x, y, subscripts, groups, ...)
>   llines(x = c(1, 6),  y = c(-24.28, 35.941667),
> lwd = 2, lty = 3, col = 4)
> }
>
> Also, it would be easier to provide tested answers if you gave a  
> reproducible example.

Here is a reproducible example in which I tried (unsuccesfully) to  
apply your suggestion:

est <- c(1:4, 3:6, 7, 9, 11, 13, 12, 15, 18, 21)
cond <- rep(c('a','b'), each = 8)
grp <- rep(c('I', 'II'), each = 4, 2)
x <- rep(c(.5, .7, .9, 1.1), 4)
upper <- est + 1
lower = est - 1
data <- data.frame(est = est, x = x, cond = cond, grp = grp, upper =  
upper, lower = lower)

xyplot(est ~ x | cond, group = grp, data = data, type = 'b',
 panel = function(x, y, subscripts, groups, ...){
 panel.superpose(c(x, rev(x)), c(upper, rev(lower)),  
subscripts, groups,
 panel.groups = 'panel.polygon', default.units =  
'native', ...)
 panel.xyplot(x, y, subscripts, groups, ...)
 }
)
# returns: "Error in validGP(list(...)) : Must specify only one of  
'font' and 'fontface'", which is incomprehensible to me.
_
Professor Michael Kubovy
University of Virginia
Department of Psychology
USPS: P.O.Box 400400Charlottesville, VA 22904-4400
Parcels:Room 102Gilmer Hall
 McCormick RoadCharlottesville, VA 22903
Office:B011+1-434-982-4729
Lab:B019+1-434-982-4751
Fax:+1-434-982-4766
WWW:http://www.people.virginia.edu/~mk9y/

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

[R] xyplot() and controlling panel.polygon()

[Michael Kubovy]

Dear R-helpers,

How do I tell panel.polygon what greoup and panel it applies to whithin
xyplot(y ~ x | c, groups = g
panel = function(x, y, groups, ...){
panel.polygon(x = xpol[c, g], y = ypol[c, g], default.units  
= 'native')
panel.xYplot(x, y, groups, ...)
llines(x = c(1, 6),  y = c(-24.28, 35.941667), lwd = 2, lty 
=  
3, col = 4)
 }
x[c, g] and y[c, g] describe the polygon I want plotted for group g  
in panel c.
_
Professor Michael Kubovy
University of Virginia
Department of Psychology
USPS: P.O.Box 400400Charlottesville, VA 22904-4400
Parcels:Room 102Gilmer Hall
 McCormick RoadCharlottesville, VA 22903
Office:B011+1-434-982-4729
Lab:B019+1-434-982-4751
Fax:+1-434-982-4766
WWW:http://www.people.virginia.edu/~mk9y/

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

[R] Bug in xYplot() with method = 'filled bands'

[Michael Kubovy]

Dear r-helpers,

I have been bitten by a cryptic comment in the help page for xYplot().

Here is some code from Madeline Bauer, one of the authors, which is  
essentailly the code on the help page for xYplot():

# This example uses the summarize function in Hmisc to
# compute the median and outer quartiles.  The outer quartiles are
# displayed using "filled bands"

set.seed(111)
dfr <- expand.grid(month=1:12, year=c(1997,1998), reps=1:100)
month <- dfr$month; year <- dfr$year
y <- abs(month-6.5) + 2*runif(length(month)) + year-1997
s <- summarize(y, llist(month,year), smedian.hilow, conf.int=.5)

# filled bands: default fill = pastel colors matching solid colors
# in superpose.line (this works differently in R)
xYplot ( Cbind ( y, Lower, Upper ) ~ month, groups=year,
  method="filled bands" , data=s, type="l")

The comment "this works differently in R" should be ammended to a  
more explicit danger sign (otherwise users like myself might find  
themselves wasting a great deal of time). There appear to be two  
differences between S+ and R here: (1) the default colors of the  
bands are not as expected (but they can be changed at will, so that's  
not a serious problem); (2) I have not been able get rid of the  
outline around the filled bands or change their colors to match the  
fill.

I would have thought that ...

xYplot ( Cbind ( y, Lower, Upper ) ~ month, groups=year,
  method="filled bands" , data=s, type="l", lwd.bands = 0,
  lty.bands = 0, col.bands = 'white')

... would do the trick.

I even tried to fix the problem in the code of xYplot:

 > panel.XYplot

function (x, y, subscripts, groups = NULL, type = if (is.function 
(method) ||
...
if (method == "filled bands") {
if (.R.)
grid.polygon(x = c(x, rev(x)), y = c(lower, rev(upper)),
gp = gpar(fill = col.fill, col = col.bands,
lty = lty.bands, lwd = lwd.bands),
default.units = "native")
else polygon(x = c(x, rev(x)), y = c(lower, rev(upper)),
col = col.fill)
}
...
}

Obviously I don't know grid well enough to fix the problem. Help  
would be much appreciated, as well as a more explicit warning on the  
xYplot() help page.

There is, of course, the possibility that there's something wrong  
with my version of R or my OS or my machine. In that case a  
replication or non-replication of this problem on a different OS  
might alert me to this.

R version 2.5.0 RC (2007-04-21 r41269)
i386-apple-darwin8.9.1

locale:
C

attached base packages:
[1] "grDevices" "datasets"  "grid"  "graphics"  "stats" "utils"
[7] "methods"   "base"

other attached packages:
 coda lme4   Matrix  car   
cacheSweave   stashR
 "0.10-7"  "0.9975-13"  "0.9975-11"  "1.2-1" 
"0.3"  "0.2-1"
 filehash   digest circular boot  gmodels  
multcompView
  "0.9-4"  "0.2.3"  "0.3-6" "1.2-27"  
"2.13.2"  "0.1-0"
  reshape  vcd   colorspace  effects 
Hmisc  acepack
  "0.7.2"  "1.0-3"   "0.95" "1.0-10"  "3.3-1" 
"1.3-2.2"
   HH multcomp  mvtnorm  JGR
iplots   JavaGD
 "1.18-3""0.991-9"  "0.7-5" "1.4-15"   
"1.0-7"  "0.3-6"
 MASS gridBase  lattice latticeExtra 
rJava   xtable
 "7.2-33"  "0.4-3" "0.15-4"  "0.1-6"  
"0.4-15"  "1.4-3"

Thanks,
Michael Kubovy

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

Re: [R] Changing appearance of band edges in xYplot()

[Michael Kubovy]

On Apr 21, 2007, at 11:08 AM, Michael Kubovy wrote:

> Dear R-helpers,
>
> I would like the bands in an xYplot() to be gray without dark edges.
> I've tried playing around with col.bands, but that seems not to
> affect the edges, although changin col.fill does change the fill as I
> would expect.
>
> I would appreciate pointers.

I forgot to give:

R version 2.5.0 RC (2007-04-19 r41248)
i386-apple-darwin8.9.1

locale:
C

attached base packages:
[1] "grDevices" "datasets"  "grid"  "graphics"  "stats" "utils"
[7] "methods"   "base"

other attached packages:
 coda  gmodelsHmisc  acepack  
lme4   Matrix
 "0.10-7" "2.13.2"  "3.3-1""1.3-2.2"  "0.9975-13"   
"0.9975-11"
cacheSweave   stashR filehash   digest  JGR
iplots
"0.3"  "0.2-1"  "0.9-4"  "0.2.3"  
"1.4-15"  "1.0-7"
   JavaGD     MASS gridBase  lattice  
latticeExtrarJava
  "0.3-6" "7.2-33"  "0.4-3" "0.15-4"  "0.1-6"  
"0.4-15"
   xtable
  "1.4-3"
_
Professor Michael Kubovy
University of Virginia
Department of Psychology
USPS: P.O.Box 400400Charlottesville, VA 22904-4400
Parcels:Room 102Gilmer Hall
 McCormick RoadCharlottesville, VA 22903
Office:B011+1-434-982-4729
Lab:B019+1-434-982-4751
Fax:+1-434-982-4766
WWW:http://www.people.virginia.edu/~mk9y/

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

[R] Changing appearance of band edges in xYplot()

[Michael Kubovy]

Dear R-helpers,

I would like the bands in an xYplot() to be gray without dark edges.  
I've tried playing around with col.bands, but that seems not to  
affect the edges, although changin col.fill does change the fill as I  
would expect.

I would appreciate pointers.

_
Professor Michael Kubovy
University of Virginia
Department of Psychology
USPS: P.O.Box 400400Charlottesville, VA 22904-4400
Parcels:Room 102Gilmer Hall
 McCormick RoadCharlottesville, VA 22903
Office:B011+1-434-982-4729
Lab:B019+1-434-982-4751
Fax:+1-434-982-4766
WWW:http://www.people.virginia.edu/~mk9y/

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

Re: [R] Positioning in xyplot

[Michael Kubovy]

Hi Deepayan,

There appears to be a small problem in the code you created to  
reposition one panel in xyplot: it plots the panels but not the  
contents. Putting together everything that you did on this problem, I  
have the following:

require(lattice)
require(grid)
## this is a safer version of current.panel.limits()
current.limits <-
function ()
{
xlim <- convertX(unit(c(0, 1), "npc"), "native", valueOnly = TRUE)
ylim <- convertY(unit(c(0, 1), "npc"), "native", valueOnly = TRUE)
if (any(!is.finite(xlim))) xlim <- c(0, 1)
if (any(!is.finite(ylim))) ylim <- c(0, 1)
list(xlim = xlim, ylim = ylim)
}
## this calls 'fun' after moving its viewport if panel.number() == 5
callAfterMoving <-
   function(fun, border = TRUE, move.x = 1, ...)
{
   if (panel.number() != 5) {  ## was == 5
   fun(...)
   if (border) grid.rect()
   }
   else {
   cpl <- current.limits()
   pushViewport(viewport(x = move.x,
 width = unit(1, "npc"),
 xscale = cpl$xlim,
 yscale = cpl$ylim,
 clip = "off"))
   fun(...)
   if (border) grid.rect()
   upViewport()
   }
}
## panel function with axes on the left:
panel.leftaxes <- function(...)
{
if (panel.number() == 5)
   panel.axis("left", outside = TRUE,
   line.col = "black")
panel.xyplot(...)
}

z <- expand.grid(x = 1:10, p = 1:5, r = 1:10)
z$y <- rnorm(nrow(z))
z$p <- factor(z$p, levels = c(1, 5, 2, 4, 3))

xyplot(y ~ x | p, z, groups = r,
   layout = c(2, 3), type = "l",
   par.settings =
   list(axis.line = list(col = "transparent"),
strip.border = list(col = "transparent")),
   panel = function(...) {
   callAfterMoving(panel.leftaxes, ...)
   },
   strip = function(...) {
   callAfterMoving(strip.default, ...)
   },
   axis = function(..., line.col, side) {
   if (side != "left" || panel.number() != 5) {
   callAfterMoving(axis.default,
   border = FALSE,
   line.col = 'black',
   side = side,
       ...)
   }
   })


I also wonder why one doesn't need to call panel.xyplot().

Thanks,
MK


_
Professor Michael Kubovy
University of Virginia
Department of Psychology
USPS: P.O.Box 400400Charlottesville, VA 22904-4400
Parcels:Room 102Gilmer Hall
 McCormick RoadCharlottesville, VA 22903
Office:B011+1-434-982-4729
Lab:B019+1-434-982-4751
Fax:+1-434-982-4766
WWW:http://www.people.virginia.edu/~mk9y/

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

Re: [R] Positioning in xyplot

[Michael Kubovy]

Wow! Thanks to Sundar and to Deepayan for this selfless effort.(I  
rather suspect they had a good time figuring this out.) The outcome  
is that I will be able to publish just the right figure for the job.  
I will be sure to aknowledge your contribtution.

Michael

On Apr 12, 2007, at 12:48 AM, Deepayan Sarkar wrote:

> On 4/11/07, Sundar Dorai-Raj <[EMAIL PROTECTED]> wrote:
>>
>> Hi, Deepayan,
>>
>> See the attached image for what your code produced. Not sure if  
>> this is
>> what you intended.
>
> Here's the correct version of callAfterMoving (I thought I had fixed
> it, but I guess I pasted the wrong thing):
>
> ## this calls 'fun' after moving its viewport if panel.number() == 5
>
> callAfterMoving <-
>   function(fun, border = TRUE, move.x = 1, ...)
> {
>   if (panel.number() != 5) {  ## was == 5
>   fun(...)
>   if (border) grid.rect()
>   }
>   else {
>   cpl <- current.limits()
>   pushViewport(viewport(x = move.x,
> width = unit(1, "npc"),
> xscale = cpl$xlim,
> yscale = cpl$ylim,
> clip = "off"))
>   fun(...)
>   if (border) grid.rect()
>   upViewport()
>   }
> }
>
> -Deepayan

__
[EMAIL PROTECTED] mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

Re: [R] Positioning in xyplot

[Michael Kubovy]

On Apr 10, 2007, at 7:07 PM, Deepayan Sarkar wrote:

> On 4/10/07, Michael Kubovy <[EMAIL PROTECTED]> wrote:
>> Dear R-helpers,
>>
>> I have an xyplot
>>
>> of the following kind:
>> xYplot(y ~ x | p, groups = factor(gg, levels = c('1', '5', '2', '4',
>> '3')),
>> strip = strip.custom(var.name = 'g', factor.levels = c(1',  
>> '5', '2',
>> '4', '3'),
>> strip.levels = T, strip.names = T, sep = ' = ',  
>> shingle.intervals =
>> NULL),
>> data = df, type = 'b', label.curves = F, layout = c(2, 3),
>> )
>>
>> Currently this puts the panels as follows
>>
>> 3
>> 2  4
>> 1  5
>>
>> I need:
>>
>>3
>> 2   4
>> 1   5
>>
>> How can I do this? Any help will be much appreciated.
>
> I don't understand what you mean. If you meant this:
>
>
>+---+
>|   |
> +---+---+
> |   |   |
> +---+---+
> |   |   |
> +---+---+
>
> then see Sundar's answer. If you meant this:
>
>  +---+
>  |   |
> +-+---+-+
> |   |   |
> +---+---+
> |   |   |
> +---+---+
>
> then you are out of luck (unless you are willing to do some low  
> level coding).

I guess I'm out of luck (I understood your two figures after turning  
them into fixed width font). At least you put me out of *this*  
misery. Thanks so much,

Michael
_
Professor Michael Kubovy
University of Virginia
Department of Psychology
USPS: P.O.Box 400400Charlottesville, VA 22904-4400
Parcels:Room 102Gilmer Hall
 McCormick RoadCharlottesville, VA 22903
Office:B011+1-434-982-4729
Lab:B019+1-434-982-4751
Fax:+1-434-982-4766
WWW:http://www.people.virginia.edu/~mk9y/



[[alternative HTML version deleted]]

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

[R] Positioning in xyplot

[Michael Kubovy]

Dear R-helpers,

I have an xyplot

of the following kind:
xYplot(y ~ x | p, groups = factor(gg, levels = c('1', '5', '2', '4',  
'3')),
strip = strip.custom(var.name = 'g', factor.levels = c(1', '5', '2',  
'4', '3'),
strip.levels = T, strip.names = T, sep = ' = ', 
shingle.intervals =  
NULL),
data = df, type = 'b', label.curves = F, layout = c(2, 3),
)

Currently this puts the panels as follows

3
2  4
1  5

I need:

   3
2   4
1   5

How can I do this? Any help will be much appreciated.

_
Professor Michael Kubovy
University of Virginia
Department of Psychology
USPS: P.O.Box 400400Charlottesville, VA 22904-4400
Parcels:Room 102Gilmer Hall
 McCormick RoadCharlottesville, VA 22903
Office:B011+1-434-982-4729
Lab:B019+1-434-982-4751
Fax:+1-434-982-4766
WWW:http://www.people.virginia.edu/~mk9y/

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

[R] Failure of mcsamp() but not mcmcsamp()

[Michael Kubovy]

Daer r-helpers,

Can anyone help with the following:

I get:

 > mcmcsamp(txt.lmer3f)
Markov Chain Monte Carlo (MCMC) output:
Start = 1
End = 1
Thinning interval = 1


But:

 > mcsamp(txt.lmer3f)
Error in as.bugs.array(sims, program = "lmer", n.iter = n.iter,  
n.burnin = n.burnin,  :
error in parameter sigma.s in parameters.to.save


I had attached the arm package:
 > library(arm)
Loading required package: R2WinBUGS

arm (Version 1.0-15, built: 2007-3-17)
Working directory is /Users/mk/Documents/projects/Mike/prediss/Exp-Text
Loading required package: foreign

 > sessionInfo()
R version 2.4.1 (2006-12-18)
i386-apple-darwin8.8.1

locale:
C

attached base packages:
[1] "tools" "grid"  "datasets"  "stats" "graphics"   
"grDevices" "utils" "methods"
[9] "base"

other attached packages:
  foreign  armR2WinBUGS coda 
Hmisc  gmodels  effects lme4
 "0.8-18" "1.0-15"  "2.0-4" "0.10-7"  "3.2-1"  
"2.13.2"  "1.0-9"  "0.9975-13"
   Matrix  car   weavercodetools
digest   xtable latticeExtra  lattice
"0.9975-11"  "1.2-1"  "1.0.1"  "0.0-3"  "0.2.3"   
"1.4-3"      "0.1-4""0.14-17"
 gridBase MASS  JGR   iplots
JavaGDrJava
  "0.4-3" "7.2-33" "1.4-15"  "1.0-5"  "0.3-6"  
"0.4-14"

Thanks,
MK
_
Professor Michael Kubovy
University of Virginia
Department of Psychology
USPS: P.O.Box 400400Charlottesville, VA 22904-4400
Parcels:Room 102Gilmer Hall
 McCormick RoadCharlottesville, VA 22903
Office:B011+1-434-982-4729
Lab:B019+1-434-982-4751
Fax:+1-434-982-4766
WWW:http://www.people.virginia.edu/~mk9y/

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

Re: [R] partial R

[Michael Kubovy]

On Apr 2, 2007, at 5:49 AM, Pedram Rowhani wrote:

> i am wondering if there is a command in R that will give me the
> partial regression coefficients

To answer your question, you could have started with
RSiteSearch('partial regression')

It's then likely that you would figured out that one way to proceed is
install.packages('car')
?cr.plots

(You may have to restart R to get the help on a newly-installed  
package.)
_________
Professor Michael Kubovy
University of Virginia
Department of Psychology
USPS: P.O.Box 400400Charlottesville, VA 22904-4400
Parcels:Room 102Gilmer Hall
 McCormick RoadCharlottesville, VA 22903
Office:B011+1-434-982-4729
Lab:B019+1-434-982-4751
Fax:+1-434-982-4766
WWW:http://www.people.virginia.edu/~mk9y/

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

[R] Why does lmList() fail when lm() doesn't?

[Michael Kubovy]

Dear r-helpers,

Can anyone suggest why lm() doesn't complain here:

summary(osss.lm1 <- lm(logOdds ~ c.setSize %in% task, data = osss))

whereas in package:nlme (and in package:lme4)

osss.lmL <- lmList(logOdds ~ c.setSize %in% task | subj, data = osss)
# Error in `contrasts<-`(`*tmp*`, value = "contr.treatment") :
#   contrasts can be applied only to factors with 2 or more levels

If it were because for each subj I have only 1 datum per cell:

with(subset(osss, osss$subj == '42'), table(setSize, task))
#   task
# setSize os ss
#  2  1  1
#  3  1  1
#  4  1  1
#  5  1  1
#  6  1  0
#  7  1  0

then

osss.lm1 <- lm(logOdds ~ c.setSize %in% task, data = subset(osss, osss 
$subj == '42'))

should fail as well, but I get

anova(osss.lm1)
# Analysis of Variance Table
#
# Response: logOdds
#Df  Sum Sq Mean Sq F value Pr(>F)
# c.setSize:task  2  6.8269  3.4134  1.9976 0.2059
# Residuals   7 11.9612  1.7087


_
Professor Michael Kubovy
University of Virginia
Department of Psychology
USPS: P.O.Box 400400Charlottesville, VA 22904-4400
Parcels:Room 102Gilmer Hall
 McCormick RoadCharlottesville, VA 22903
Office:B011+1-434-982-4729
Lab:B019+1-434-982-4751
Fax:+1-434-982-4766
WWW:http://www.people.virginia.edu/~mk9y/

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

Re: [R] Doing partial-f test for stepwise regression

[Michael Kubovy]

On Apr 1, 2007, at 1:54 AM, [EMAIL PROTECTED] wrote:

> Hello all,
> I am trying to figure out an optimal linear model by using stepwise
> regression which requires partial f-test, I did some Googling on the
> Internet and realised that someone seemed to ask the question before:
>
> Jim Milks <[EMAIL PROTECTED]> writes:
>> Dear all:
>>
>> I have a regression model that has collinearity problems (between
>> three regressor variables). I need a F-test that will allow me to
>> compare between full (with all variables) and partial models (minus
>> 1=< variables). The general F-test formula I'm using is:
>>
>> F = {[SS(full model) - SS(reduced model)] / (#variables taken out)} /
>> MSS(full model)
>>
>> Unfortunately, the ANOVA table parses the SS and MSS between the
>> variables and does not give the statistics for the regression  
>> model as
>> a whole, otherwise I'd do this by hand.
>>
>> So, really, I have two questions: 1) Can I just add up all the SS and
>> MSS for all the variables to get the model SS and MSS and 2) Are
>> there any functions or packages I can use to calculate the F- 
>> statistic?
>> Just use anova(model1, model2).
>> (One potential catch: Make sure that both models are fitted to the  
>> same
>> data set. Missing values in predictors may interfere.)
>
> However, in the answer provided by Mr. Peter Dalgaard,(use
> anova(model1,model2) I could not understand what model1 and model2 are
> supposed to referring to, which one is supposedly to be the full  
> model and
> which one is to be the partial model? Or it does not matter?

You can tell which is which by looking at the degrees of freedom.

_
Professor Michael Kubovy
University of Virginia
Department of Psychology
USPS: P.O.Box 400400Charlottesville, VA 22904-4400
Parcels:Room 102Gilmer Hall
 McCormick RoadCharlottesville, VA 22903
Office:B011+1-434-982-4729
Lab:B019+1-434-982-4751
Fax:+1-434-982-4766
WWW:http://www.people.virginia.edu/~mk9y/

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

Re: [R] Newbie: Combn and scripting

[Michael Kubovy]

x <- matrix(c('a', 'b', 'c', 'd', 1:4, 'g', 3, 6, 't'), nrow = 3,  
byrow = T)
comb <- vector('list', 3)
for(i in 1:3) comb[[i]] <- combn(x[i,], 3)

On Mar 27, 2007, at 1:55 AM, [EMAIL PROTECTED] wrote:

> I have just installed my R 2.4 (windows) as a test trying to load a  
> data
> frame and run combn() for each line into another file. How do I do  
> this?
>
> data.csv:
> a,b,c,d
> 1,2,3.4
> g,3,6,t
> etc
>
> x=data.csv, m=3



_
Professor Michael Kubovy
University of Virginia
Department of Psychology
USPS: P.O.Box 400400Charlottesville, VA 22904-4400
Parcels:Room 102Gilmer Hall
 McCormick RoadCharlottesville, VA 22903
Office:B011+1-434-982-4729
Lab:B019+1-434-982-4751
Fax:+1-434-982-4766
WWW:http://www.people.virginia.edu/~mk9y/

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

Re: [R] distribution graph

[Michael Kubovy]

?violinplot (You need to install the UsingR package first.)

On Mar 23, 2007, at 4:06 AM, Plessen, Christian von wrote:

> I am looking for a way to produce a "distribution graph" as in the  
> example:
>
> (http://cecsweb.dartmouth.edu/release1.1/datatools/dgraph.php? 
> year=2003&geotype=STD_HRR&event=A01_DIS&eventtype=UTIL
>
> Anybody who can help?

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

Re: [R] package:AlgDesign and .Random.seed

[Michael Kubovy]

On Mar 21, 2007, at 4:16 AM, Uwe Ligges wrote:

> Michael Kubovy wrote:
>> Dear r-helpers,
>> Could you please help me solve the following problem: When I run
>> require(AlgDesign)
>> trt <- LETTERS[1:5]
>> blk <- 10
>> trtblk <- 3
>> BIB <- optBlock(~., withinData = trt, blocksizes = rep(trtblk, blk))
>> In response to the last command, R complains:
>> Error in optBlock(~., withinData = trt, blocksizes = rep(trtblk,  
>> blk)) :
>>  object ".Random.seed" not found
>> The documentation of optBlock() in AlgDesign doesn't say that I   
>> needed to set .Random.seed. I thought it was initiated  
>> automatically  at the beginning of a session. What am I missing?
>
>
> The first line in that function is
> seed <- .Random.seed
> but .Random.seed is generated at the first use of R's RNG, hence  
> maybe later. This means the function contains a bug which you  
> should report to the package maintainer, please.
>
> Best,
> Uwe Ligges

Bob Wheeler's response:

> From: Bob Wheeler <[EMAIL PROTECTED]>
> Date: March 21, 2007 9:19:29 AM EDT
> To: Michael Kubovy <[EMAIL PROTECTED]>
> Subject: Re:
>
> Each workspace in R requires you to set a random seed to start. You  
> have not done this. It is an R artifact, and has nothing to do with  
> AlgDesign.
> -- 
> Bob Wheeler --- http://www.bobwheeler.com/
>ECHIP, Inc. --- Randomness comes in bunches.

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

[R] package:AlgDesign and .Random.seed

[Michael Kubovy]

Dear r-helpers,

Could you please help me solve the following problem: When I run

require(AlgDesign)
trt <- LETTERS[1:5]
blk <- 10
trtblk <- 3
BIB <- optBlock(~., withinData = trt, blocksizes = rep(trtblk, blk))

In response to the last command, R complains:

Error in optBlock(~., withinData = trt, blocksizes = rep(trtblk, blk)) :
object ".Random.seed" not found

The documentation of optBlock() in AlgDesign doesn't say that I  
needed to set .Random.seed. I thought it was initiated automatically  
at the beginning of a session. What am I missing?

 > sessionInfo()
R version 2.4.1 (2006-12-18)
i386-apple-darwin8.8.1

locale:
C

attached base packages:
[1] "grid"  "datasets"  "stats" "graphics"  "grDevices"  
"utils" "methods"   "base"

other attached packages:
AlgDesign   xtable latticeExtra  lattice  
gridBase MASS  JGR   iplots
  "1.0-7"  "1.4-3"  "0.1-4""0.14-16"  "0.4-3"  
"7.2-32" "1.4-15"  "1.0-5"
   JavaGDrJava
  "0.3-6" "0.4-14"





_
Professor Michael Kubovy
University of Virginia
Department of Psychology
USPS: P.O.Box 400400Charlottesville, VA 22904-4400
Parcels:Room 102Gilmer Hall
 McCormick RoadCharlottesville, VA 22903
Office:B011+1-434-982-4729
Lab:B019+1-434-982-4751
Fax:+1-434-982-4766
WWW:http://www.people.virginia.edu/~mk9y/

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

[R] Help with xYplot() in package Hmisc

[Michael Kubovy]

Dear r-helpers,

I would like to add three fitted lines to this plot. I don't know how  
to include three panel.abline() to the call to xYplot. I would like  
the colors of the lines correspond to the color of the dots and CI bars.

x <- c(1:4, 1:4, 3:6)
y <- c(4:1, 6:3, 5:2)
e <- runif(12) - .5
y <- y + e
ll <- y - 1
ul <- y + 1
g <- c(rep(1, 4), rep(2, 4), rep(3, 4))
df <- data.frame(x, y, ll, ul, g)
xYplot(Cbind(y, ll, ul) ~ x, groups = g, type = 'p', data = df)
_________
Professor Michael Kubovy
University of Virginia
Department of Psychology
USPS: P.O.Box 400400Charlottesville, VA 22904-4400
Parcels:Room 102Gilmer Hall
 McCormick RoadCharlottesville, VA 22903
Office:B011+1-434-982-4729
Lab:B019+1-434-982-4751
Fax:+1-434-982-4766
WWW:http://www.people.virginia.edu/~mk9y/

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

Re: [R] logging mouse clicks

[Michael Kubovy]

Hi Seth,

On Mar 10, 2007, at 11:55 PM, Seth Roberts wrote:

> How can I use R to record the time of a mouse click?

Perhaps ievent.wait() in iplots? (I haven't used this.) I found all  
sorts of functions (most of which I couldn't understand) by ?event.
_____
Professor Michael Kubovy
University of Virginia
Department of Psychology
USPS: P.O.Box 400400Charlottesville, VA 22904-4400
Parcels:Room 102Gilmer Hall
 McCormick RoadCharlottesville, VA 22903
Office:B011+1-434-982-4729
Lab:B019+1-434-982-4751
Fax:+1-434-982-4766
WWW:http://www.people.virginia.edu/~mk9y/

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

Re: [R] Table Construction from calculations

[Michael Kubovy]

On Mar 9, 2007, at 9:23 PM, Seth Imhoff wrote:

> I am trying to create a table of values by adding  pairs of  
> vectors, but
> am running into some problems.  The problem is best expressed by a
> simple example.
>
> Starting with a data table "basis":
>   atom   x   y   z
> 1   Cu 0.0 0.0 0.0
> 2   Cu 0.5 0.5 0.5
>
> I want to add 0.5 0.5 0.5  (and also the 0 0 0 but it wouldn't change
> the values below so I won't refer to it in the rest of the example)  
> to a
> list of vectors in the form of:
>> latpoints
>V1 V2 V3
> 1   0  0  0
> 2   0  0  1
> 3   0  0  2
> 4   0  0  3
> 5   0  1  1
>
> so that I end up with a table such as:
> V1  V2  V3
> 0.5 0.5 0.5
> 0.5 0.5 1.5
> 0.5 0.5 2.5
> 0.5 0.5 3.5
> 0.5 1.5 1.5
>
> I've tried many variations on the following: (not just cat, but  
> most of
> the data/data.table options)
>
>  test = for(i in 1:5) {cat(basis[1,2:4] + latticemultipliers[i,],
> append=TRUE)}
>
> However, I either end up with an error telling me that cat doesn't
> handle "type 'list' " or with a table with length of 1 such as:
>  xyz
> 2  0.5 1.5 1.5

Is this what you want?

 > (latpoints <- data.frame(matrix(c(0,0,0,0,0,1,0,0,2,0,0,3,0,1,1),  
nrow = 5, byrow = T)))
   X1 X2 X3
1  0  0  0
2  0  0  1
3  0  0  2
4  0  0  3
5  0  1  1
 > (latpoints <- latpoints + c(0.5, 0.5, 0.5))
X1  X2  X3
1 0.5 0.5 0.5
2 0.5 0.5 1.5
3 0.5 0.5 2.5
4 0.5 0.5 3.5
5 0.5 1.5 1.5

This is an important feature of R called "vectorization" (see, .e.g,  
cran.r-project.org/doc/contrib/Paradis-rdebuts_en.pdf or  
www.ms.washington.edu/stat390/winter07/R_primer.pdf) which allows you  
do avoid writing loops.
_
Professor Michael Kubovy
University of Virginia
Department of Psychology
USPS: P.O.Box 400400Charlottesville, VA 22904-4400
Parcels:Room 102Gilmer Hall
 McCormick RoadCharlottesville, VA 22903
Office:B011+1-434-982-4729
Lab:B019+1-434-982-4751
Fax:+1-434-982-4766
WWW:http://www.people.virginia.edu/~mk9y/

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

Re: [R] Extracting the p of F statistics from lm

[Michael Kubovy]

On Mar 9, 2007, at 11:18 AM, Cressoni, Massimo ((NIH/NHLBI)) [F] wrote:

> I need to extract the p value from a ANOVA done with lm model
>
> fitting <- lm(var ~ group)
> Sfitting <- summary(fitting)
>
> Sfitting[10][1] gives the F value and the degrees of freedom but I  
> am not able to get the
> p value.

try
Sfitting[4]$coefficients[,4]

I'm not sure that this is the best way, but it works with the example  
for lm()
 > summary(lm.D9)[4]$coefficients[,4]
# (Intercept) groupTrt
# 9.547128e-15 2.490232e-01

_________
Professor Michael Kubovy
University of Virginia
Department of Psychology
USPS: P.O.Box 400400Charlottesville, VA 22904-4400
Parcels:Room 102Gilmer Hall
 McCormick RoadCharlottesville, VA 22903
Office:B011+1-434-982-4729
Lab:B019+1-434-982-4751
Fax:+1-434-982-4766
WWW:http://www.people.virginia.edu/~mk9y/

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

Re: [R] Is the gmodels package being maintained?

[Michael Kubovy]

Hi,

Finding his email address was not immediate, but I finally did, and  
did bring the problem to Greg's attention @ rochester, and the  
message didn't bounce this time.

On Mar 9, 2007, at 7:43 AM, Peter Dalgaard wrote:

> Michael Kubovy wrote:
>> Dear r-helpers,
>>
>> I sent  a cc of a recent message about a problem with ci.lmer() in
>> the gmodels package to the author (Gregory R Warnes), and the message
>> bounced. If the author or someone else is maintaining this package or
>> this function, would you kindly supplement the author's name and/or
>> address with a current maintainer and/or provide a current email
>> address?
>>
> Haven't heard that Greg should be out of circulation. You might try  
> the
> address from his homepage:
>
> [EMAIL PROTECTED]
_
Professor Michael Kubovy
University of Virginia
Department of Psychology
USPS: P.O.Box 400400Charlottesville, VA 22904-4400
Parcels:Room 102Gilmer Hall
 McCormick RoadCharlottesville, VA 22903
Office:B011+1-434-982-4729
Lab:B019+1-434-982-4751
Fax:+1-434-982-4766
WWW:http://www.people.virginia.edu/~mk9y/

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

[R] Is the gmodels package being maintained?

[Michael Kubovy]

Dear r-helpers,

I sent  a cc of a recent message about a problem with ci.lmer() in  
the gmodels package to the author (Gregory R Warnes), and the message  
bounced. If the author or someone else is maintaining this package or  
this function, would you kindly supplement the author's name and/or  
address with a current maintainer and/or provide a current email  
address?
_
Professor Michael Kubovy
University of Virginia
Department of Psychology
USPS: P.O.Box 400400Charlottesville, VA 22904-4400
Parcels:Room 102Gilmer Hall
 McCormick RoadCharlottesville, VA 22903
Office:B011+1-434-982-4729
Lab:B019+1-434-982-4751
Fax:+1-434-982-4766
WWW:http://www.people.virginia.edu/~mk9y/

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

[R] Problem with ci.lmer() in package:gmodels

[Michael Kubovy]

Dear Friends,

Please note that in the following CI lower > CI higher:

 > require(lmer)
 > require(gmodels)
 > fm2 <- lmer(Reaction ~ Days + (1|Subject) + (0+Days|Subject),  
sleepstudy)
 > ci(fm2)
  Estimate  CI lower   CI upper Std. Error p-value
(Intercept) 251.66693 266.06895 238.630280   7.056447   0
Days 10.52773  13.63372   7.389946   1.646900   0



_____
Professor Michael Kubovy
University of Virginia
Department of Psychology
USPS: P.O.Box 400400Charlottesville, VA 22904-4400
Parcels:Room 102Gilmer Hall
 McCormick RoadCharlottesville, VA 22903
Office:B011+1-434-982-4729
Lab:B019+1-434-982-4751
Fax:+1-434-982-4766
WWW:http://www.people.virginia.edu/~mk9y/

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

Re: [R] Failure to run mcsamp() in package arm

[Michael Kubovy]

Andrew Robinson has gently chided me for not including more  
information. So here goes:
R version 2.4.1 (2006-12-18)
powerpc-apple-darwin8.8.0

locale:
C

attached base packages:
[1] "stats" "graphics"  "grDevices" "utils" "methods"   "base"

other attached packages:
 foreign car arm   R2WinBUGSlme4   
Matrix lattice
"0.8-18" "1.2-1""1.0-13" "2.0-4" "0.9975-13"  
"0.9975-11"   "0.14-16"
MASS JGR  iplots  JavaGD   rJava
"7.2-32""1.4-15" "1.0-5" "0.3-5""0.4-14"

On Mar 7, 2007, at 4:30 PM, Michael Kubovy wrote:

> More problems. If I run
> sim(fm1 <- lmer(Reaction ~ Days + (Days|Subject), sleepstudy))
> from the lmer() help page.
>
> I get the error
> Error in mvrnorm(n.sims, bhat[j, ], V.beta) :
>   'Sigma' is not positive definite
>
> On Mar 7, 2007, at 1:30 PM, Michael Kubovy wrote:
>
>> Dear r-helpers,
>>
>> I can run the examples on the mcsamp help page. For example:
>> 
>>> M1 <- lmer (y1 ~ x + (1|group))
>>> (M1.sim <- mcsamp (M1))
>> fit using lmer,
>> 3 chains, each with 1000 iterations (first 500 discarded)
>> n.sims = 1500 iterations saved
>>mean  sd 2.5%  25%  50%  75% 97.5% Rhat
>> n.eff
>> beta.(Intercept)   0.1 0.7 -1.2 -0.3  0.1  0.5   1.4  1.0
>> 1500
>> beta.x 2.5 0.4  1.7  2.2  2.5  2.7   3.2  1.0
>> 1500
>> sigma.y3.8 0.3  3.3  3.6  3.7  3.9   4.3
>> 1.061
>> sigma.grop.(In)1.5 0.8  0.0  1.0  1.4  1.9   3.3
>> 1.412
>> eta.group.(Intercept)[1]   0.0 1.0 -2.1 -0.5  0.0  0.6   2.0  1.0
>> 1500
>> eta.group.(Intercept)[2]   1.0 1.1 -0.9  0.2  0.9  1.7   3.4
>> 1.059
>> eta.group.(Intercept)[3]  -1.3 1.2 -4.0 -2.0 -1.3 -0.4   0.5
>> 1.066
>> eta.group.(Intercept)[4]   1.3 1.1 -0.6  0.4  1.1  2.0   3.7
>> 1.143
>> eta.group.(Intercept)[5]  -0.7 1.0 -3.0 -1.4 -0.6  0.0   1.2  1.0
>> 120
>> eta.group.(Intercept)[6]   1.5 1.2 -0.3  0.6  1.4  2.2   4.0
>> 1.049
>> eta.group.(Intercept)[7]   0.3 1.0 -1.7 -0.3  0.1  0.8   2.5  1.0
>> 440
>> eta.group.(Intercept)[8]  -1.6 1.2 -4.0 -2.4 -1.5 -0.6   0.3
>> 1.141
>> eta.group.(Intercept)[9]   0.4 1.0 -1.6 -0.2  0.2  0.9   2.7  1.0
>> 180
>> eta.group.(Intercept)[10] -1.0 1.1 -3.3 -1.6 -0.9 -0.2   0.8
>> 1.086
>>
>> For each parameter, n.eff is a crude measure of effective sample  
>> size,
>> and Rhat is the potential scale reduction factor (at convergence,
>> Rhat=1).
>> 
>> But when I try to do this with my own data I get an error:
>> 
>>> display(e7.lmer2)
>> lmer(formula = baLO ~ I(baRatio - 0.985) + delta + (1 + I(baRatio -
>> 0.985) + delta | subject), data = e7)
>> coef.est coef.se
>> (Intercept)-0.19 0.06
>> I(baRatio - 0.985) -4.95 0.74
>> delta   0.41 0.06
>> Error terms:
>> Groups   Name   Std.Dev. Corr
>> subject  (Intercept)0.13
>>I(baRatio - 0.985) 2.57  0.45
>>delta  0.22 -0.12 -0.94
>> Residual0.39
>> number of obs: 494, groups: subject, 13
>> deviance = 551.4
>>
>>> e7.sim <- mcsamp(e7.lmer2)
>> Error in as.bugs.array(sims, program = "lmer", n.iter = n.iter,
>> n.burnin = n.burnin,  :
>>  error in parameter sigma. in parameters.to.save
>> 
>> I would appreciate a pointer to what the problem might be.
>> _
>> Professor Michael Kubovy
>> University of Virginia
>> Department of Psychology
>> USPS: P.O.Box 400400Charlottesville, VA 22904-4400
>> Parcels:Room 102Gilmer Hall
>>  McCormick RoadCharlottesville, VA 22903
>> Office:B011+1-434-982-4729
>> Lab:B019+1-434-982-4751
>> Fax:+1-434-982-4766
>> WWW:http://www.people.virginia.edu/~mk9y/
>>
>> __
>> R-help@stat.math.ethz.ch mailing list
>> https://stat.ethz.ch/mailman/listinfo/r-help
>> PLEASE do read the posting guide http://www.R-project.org/posting-
>> guide.html
>> and provide commented, minimal, self-contained, reproducible code.
>
> __
> R-help@stat.math.ethz.ch mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide http://www.R-project.org/posting- 
> guide.html
> and provide commented, minimal, self-contained, reproducible code.

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

Re: [R] Failure to run mcsamp() in package arm

[Michael Kubovy]

More problems. If I run
sim(fm1 <- lmer(Reaction ~ Days + (Days|Subject), sleepstudy))
from the lmer() help page.

I get the error
Error in mvrnorm(n.sims, bhat[j, ], V.beta) :
'Sigma' is not positive definite

On Mar 7, 2007, at 1:30 PM, Michael Kubovy wrote:

> Dear r-helpers,
>
> I can run the examples on the mcsamp help page. For example:
> 
>> M1 <- lmer (y1 ~ x + (1|group))
>> (M1.sim <- mcsamp (M1))
> fit using lmer,
> 3 chains, each with 1000 iterations (first 500 discarded)
> n.sims = 1500 iterations saved
>mean  sd 2.5%  25%  50%  75% 97.5% Rhat  
> n.eff
> beta.(Intercept)   0.1 0.7 -1.2 -0.3  0.1  0.5   1.4  1.0   
> 1500
> beta.x 2.5 0.4  1.7  2.2  2.5  2.7   3.2  1.0   
> 1500
> sigma.y3.8 0.3  3.3  3.6  3.7  3.9   4.3   
> 1.061
> sigma.grop.(In)1.5 0.8  0.0  1.0  1.4  1.9   3.3   
> 1.412
> eta.group.(Intercept)[1]   0.0 1.0 -2.1 -0.5  0.0  0.6   2.0  1.0   
> 1500
> eta.group.(Intercept)[2]   1.0 1.1 -0.9  0.2  0.9  1.7   3.4   
> 1.059
> eta.group.(Intercept)[3]  -1.3 1.2 -4.0 -2.0 -1.3 -0.4   0.5   
> 1.066
> eta.group.(Intercept)[4]   1.3 1.1 -0.6  0.4  1.1  2.0   3.7   
> 1.143
> eta.group.(Intercept)[5]  -0.7 1.0 -3.0 -1.4 -0.6  0.0   1.2  1.0
> 120
> eta.group.(Intercept)[6]   1.5 1.2 -0.3  0.6  1.4  2.2   4.0   
> 1.049
> eta.group.(Intercept)[7]   0.3 1.0 -1.7 -0.3  0.1  0.8   2.5  1.0
> 440
> eta.group.(Intercept)[8]  -1.6 1.2 -4.0 -2.4 -1.5 -0.6   0.3   
> 1.141
> eta.group.(Intercept)[9]   0.4 1.0 -1.6 -0.2  0.2  0.9   2.7  1.0
> 180
> eta.group.(Intercept)[10] -1.0 1.1 -3.3 -1.6 -0.9 -0.2   0.8   
> 1.086
>
> For each parameter, n.eff is a crude measure of effective sample size,
> and Rhat is the potential scale reduction factor (at convergence,
> Rhat=1).
> 
> But when I try to do this with my own data I get an error:
> 
>> display(e7.lmer2)
> lmer(formula = baLO ~ I(baRatio - 0.985) + delta + (1 + I(baRatio -
> 0.985) + delta | subject), data = e7)
> coef.est coef.se
> (Intercept)-0.19 0.06
> I(baRatio - 0.985) -4.95 0.74
> delta   0.41 0.06
> Error terms:
> Groups   Name   Std.Dev. Corr
> subject  (Intercept)0.13
>I(baRatio - 0.985) 2.57  0.45
>delta  0.22 -0.12 -0.94
> Residual0.39
> number of obs: 494, groups: subject, 13
> deviance = 551.4
>
>> e7.sim <- mcsamp(e7.lmer2)
> Error in as.bugs.array(sims, program = "lmer", n.iter = n.iter,
> n.burnin = n.burnin,  :
>   error in parameter sigma. in parameters.to.save
> 
> I would appreciate a pointer to what the problem might be.
> _
> Professor Michael Kubovy
> University of Virginia
> Department of Psychology
> USPS: P.O.Box 400400Charlottesville, VA 22904-4400
> Parcels:Room 102Gilmer Hall
>  McCormick RoadCharlottesville, VA 22903
> Office:B011+1-434-982-4729
> Lab:B019+1-434-982-4751
> Fax:+1-434-982-4766
> WWW:http://www.people.virginia.edu/~mk9y/
>
> __
> R-help@stat.math.ethz.ch mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide http://www.R-project.org/posting- 
> guide.html
> and provide commented, minimal, self-contained, reproducible code.

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

[R] Failure to run mcsamp() in package arm

[Michael Kubovy]

Dear r-helpers,

I can run the examples on the mcsamp help page. For example:

 > M1 <- lmer (y1 ~ x + (1|group))
 > (M1.sim <- mcsamp (M1))
fit using lmer,
3 chains, each with 1000 iterations (first 500 discarded)
n.sims = 1500 iterations saved
   mean  sd 2.5%  25%  50%  75% 97.5% Rhat n.eff
beta.(Intercept)   0.1 0.7 -1.2 -0.3  0.1  0.5   1.4  1.0  1500
beta.x 2.5 0.4  1.7  2.2  2.5  2.7   3.2  1.0  1500
sigma.y3.8 0.3  3.3  3.6  3.7  3.9   4.3  1.061
sigma.grop.(In)1.5 0.8  0.0  1.0  1.4  1.9   3.3  1.412
eta.group.(Intercept)[1]   0.0 1.0 -2.1 -0.5  0.0  0.6   2.0  1.0  1500
eta.group.(Intercept)[2]   1.0 1.1 -0.9  0.2  0.9  1.7   3.4  1.059
eta.group.(Intercept)[3]  -1.3 1.2 -4.0 -2.0 -1.3 -0.4   0.5  1.066
eta.group.(Intercept)[4]   1.3 1.1 -0.6  0.4  1.1  2.0   3.7  1.143
eta.group.(Intercept)[5]  -0.7 1.0 -3.0 -1.4 -0.6  0.0   1.2  1.0   120
eta.group.(Intercept)[6]   1.5 1.2 -0.3  0.6  1.4  2.2   4.0  1.049
eta.group.(Intercept)[7]   0.3 1.0 -1.7 -0.3  0.1  0.8   2.5  1.0   440
eta.group.(Intercept)[8]  -1.6 1.2 -4.0 -2.4 -1.5 -0.6   0.3  1.141
eta.group.(Intercept)[9]   0.4 1.0 -1.6 -0.2  0.2  0.9   2.7  1.0   180
eta.group.(Intercept)[10] -1.0 1.1 -3.3 -1.6 -0.9 -0.2   0.8  1.086

For each parameter, n.eff is a crude measure of effective sample size,
and Rhat is the potential scale reduction factor (at convergence,  
Rhat=1).

But when I try to do this with my own data I get an error:

 > display(e7.lmer2)
lmer(formula = baLO ~ I(baRatio - 0.985) + delta + (1 + I(baRatio -   
0.985) + delta | subject), data = e7)
coef.est coef.se
(Intercept)-0.19 0.06
I(baRatio - 0.985) -4.95 0.74
delta   0.41 0.06
Error terms:
Groups   Name   Std.Dev. Corr
subject  (Intercept)0.13
   I(baRatio - 0.985) 2.57  0.45
   delta  0.22 -0.12 -0.94
Residual0.39
number of obs: 494, groups: subject, 13
deviance = 551.4

 > e7.sim <- mcsamp(e7.lmer2)
Error in as.bugs.array(sims, program = "lmer", n.iter = n.iter,  
n.burnin = n.burnin,  :
error in parameter sigma. in parameters.to.save

I would appreciate a pointer to what the problem might be.
_________
Professor Michael Kubovy
University of Virginia
Department of Psychology
USPS: P.O.Box 400400Charlottesville, VA 22904-4400
Parcels:Room 102Gilmer Hall
 McCormick RoadCharlottesville, VA 22903
Office:B011+1-434-982-4729
Lab:B019+1-434-982-4751
Fax:+1-434-982-4766
WWW:http://www.people.virginia.edu/~mk9y/

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

Re: [R] Question

[Michael Kubovy]

> How can I evaluate two or more expressions to return two or more  
> columns?
>
> eval(expression(with(bd,Var1*100),with(bd,Var2*200)))
>
> The execution is always for the last expression. I can't execute  
> the two expressions at the same time. It is possible?

How about something like this

with(bd, c(Var1 * 100, Var2 * 200))
_________
Professor Michael Kubovy
University of Virginia
Department of Psychology
USPS: P.O.Box 400400Charlottesville, VA 22904-4400
Parcels:Room 102Gilmer Hall
 McCormick RoadCharlottesville, VA 22903
Office:B011+1-434-982-4729
Lab:B019+1-434-982-4751
Fax:+1-434-982-4766
WWW:http://www.people.virginia.edu/~mk9y/

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

Re: [R] How to override ordering of panels in xyplot()

[Michael Kubovy]

On Mar 3, 2007, at 3:19 PM, Deepayan Sarkar wrote:

> On 3/3/07, Michael Kubovy <[EMAIL PROTECTED]> wrote:
>> Dear r-helpers,
>>
>> I'm conditioning an xyplot on a variable whose levels are'low',  
>> 'med', 'high'. How do I override the alphabetical ordering for the  
>> panels of the plot?
>
> This has less to do with xyplot and more to do with the default of  
> the 'levels' argument in the factor() function. Just make sure the  
> levels are in the right order in your data when xyplot is called.

Unless one makes the factor ordered, reordering the levels of a  
factor does not seem to be a trivial matter. The only R function I've  
found that makes it easy is reorder_factor() (package:reshape). Or am  
I missing something?
_
Professor Michael Kubovy
University of Virginia
Department of Psychology
USPS: P.O.Box 400400Charlottesville, VA 22904-4400
Parcels:Room 102Gilmer Hall
 McCormick RoadCharlottesville, VA 22903
Office:B011+1-434-982-4729
Lab:B019+1-434-982-4751
Fax:+1-434-982-4766
WWW:http://www.people.virginia.edu/~mk9y/

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

[R] How to override ordering of panels in xyplot()

[Michael Kubovy]

Dear r-helpers,

I'm conditioning an xyplot on a variable whose levels are'low',  
'med', 'high'. How do I override the alphabetical ordering for the  
panels of the plot?
_________
Professor Michael Kubovy
University of Virginia
Department of Psychology
USPS: P.O.Box 400400Charlottesville, VA 22904-4400
Parcels:Room 102Gilmer Hall
 McCormick RoadCharlottesville, VA 22903
Office:B011+1-434-982-4729
Lab:B019+1-434-982-4751
Fax:+1-434-982-4766
WWW:http://www.people.virginia.edu/~mk9y/

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

[R] Help with paste()

[Michael Kubovy]

Dear r-helpers,

Could you please tell me what's missing:
rbind(paste('txt.est',1:24, sep = ''))
txt.est1, ... txt.est24 are vectors that I wish to rbind.
_________
Professor Michael Kubovy
University of Virginia
Department of Psychology
USPS: P.O.Box 400400Charlottesville, VA 22904-4400
Parcels:Room 102Gilmer Hall
 McCormick RoadCharlottesville, VA 22903
Office:B011+1-434-982-4729
Lab:B019+1-434-982-4751
Fax:+1-434-982-4766
WWW:http://www.people.virginia.edu/~mk9y/

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

Re: [R] mtext bold font problem

[Michael Kubovy]

On Feb 24, 2007, at 6:12 PM, Joseph Retzer wrote:

> I'd like to create an mtext title which has certain words bold, the  
> rest not bold.

x <- 1:10
y <- x^2
plot(x, y, main = expression(paste(
 "Overall Satisfaction Call",
 bold(" Difference "),
 "vs. Overall Satisfaction Rep.",
 bold(" Difference")))
)

_
Professor Michael Kubovy
University of Virginia
Department of Psychology
USPS: P.O.Box 400400Charlottesville, VA 22904-4400
Parcels:Room 102Gilmer Hall
 McCormick RoadCharlottesville, VA 22903
Office:B011+1-434-982-4729
Lab:B019+1-434-982-4751
Fax:+1-434-982-4766
WWW:http://www.people.virginia.edu/~mk9y/

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

Re: [R] Wrinting integers in a matrix faile

[Michael Kubovy]

On Feb 22, 2007, at 8:36 PM, Mohsen Jafarikia wrote:

> I am using the following program to get the p-value of some numbers
> (column 'LR' of the data.dat file). I want to write the 1st and 2nd
> column of the output file (data.out) as an integer while the program
> change them to double. Could anybody please tell me how I can write
> the code which writes the values of the first two columns as integer?
>
> library ('MASS')
> MP<-read.table(file='data.dat')
> names(MP)<-c('B','R','S','L','LR','Q')
> a<-as.matrix((1-pchisq(MP$LR, df=1)))
> b<-cbind(MP$B,MP$R,a,MP$S,MP$L,MP$LR,MP$Q)
> write.matrix(b,'data.out')

?Round
_
Professor Michael Kubovy
University of Virginia
Department of Psychology
USPS: P.O.Box 400400Charlottesville, VA 22904-4400
Parcels:Room 102Gilmer Hall
 McCormick RoadCharlottesville, VA 22903
Office:B011+1-434-982-4729
Lab:B019+1-434-982-4751
Fax:+1-434-982-4766
WWW:http://www.people.virginia.edu/~mk9y/

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

Re: [R] question about boxplot

[Michael Kubovy]

On Feb 22, 2007, at 2:56 PM, Smith, Phil ((CDC/CCID/NCIRD)) wrote:

> boxplot( p.prop ~ R  + bins )
>
> This command makes nice boxplot for the factor "R" crossed with the  
> factor "bins".
>
> I am having alot of trouble getting control of the labels on the X  
> axis. I want to control it more by specifying what the labels are,  
> controling the 'size' of those labels (by using cex), and then  
> control the rotation of the character strings of those labels (by  
> using srt or crt).
>
> There is a "names" argument to boxplot, but I haven't had much luck  
> controlling what it prints, the size of the font, and the character  
> rotation.

Would you consider an easy way out---an alternative with reasonable  
defaults?

data(ToothGrowth)
bwplot(dose ~ len | supp, ToothGrowth)
_
Professor Michael Kubovy
University of Virginia
Department of Psychology
USPS: P.O.Box 400400Charlottesville, VA 22904-4400
Parcels:Room 102Gilmer Hall
 McCormick RoadCharlottesville, VA 22903
Office:B011+1-434-982-4729
Lab:B019+1-434-982-4751
Fax:+1-434-982-4766
WWW:http://www.people.virginia.edu/~mk9y/

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

[R] Help with xlab and ylab distance from axes

[Michael Kubovy]

Dear r-helpers,

In basic graphics, I have a figure with x and y axes suppresed. I  
would like to move the xlab and the ylab closer to the axes. Do I  
have to use mtext()?
_
Professor Michael Kubovy
University of Virginia
Department of Psychology
USPS: P.O.Box 400400Charlottesville, VA 22904-4400
Parcels:Room 102Gilmer Hall
 McCormick RoadCharlottesville, VA 22903
Office:B011+1-434-982-4729
Lab:B019+1-434-982-4751
Fax:+1-434-982-4766
WWW:http://www.people.virginia.edu/~mk9y/

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

Re: [R] Suppressing \newcolumntype{} declaration in latex.default() in Hmisc

[Michael Kubovy]

Dear r-helpers,

I didn't include a clear question in my previous posting. Here is a  
better version:

When I use latex() on a matrix, I set the option dcolumn = T. As a   
result, in the *.tex file before each tabular I get
\newcolumntype{.}{D{.}{.}{-1}}
The LaTeX compiler complains about these multiple redeclarations. Is  
there a way to suppress this?
_
Professor Michael Kubovy
University of Virginia
Department of Psychology
USPS: P.O.Box 400400Charlottesville, VA 22904-4400
Parcels:Room 102Gilmer Hall
 McCormick RoadCharlottesville, VA 22903
Office:B011+1-434-982-4729
Lab:B019+1-434-982-4751
Fax:+1-434-982-4766
WWW:http://www.people.virginia.edu/~mk9y/

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

[R] Suppressing \newcolumntype{} declaration in latex.default() in Hmisc

[Michael Kubovy]

Dear r-helpers,

When use latex() on a matrix, I set the option dcolumn = T. As a  
result, in front of each tabular I get
\newcolumntype{.}{D{.}{.}{-1}}
The LaTeX compiler complains about these multiple redeclarations.
_
Professor Michael Kubovy
University of Virginia
Department of Psychology
USPS: P.O.Box 400400Charlottesville, VA 22904-4400
Parcels:Room 102Gilmer Hall
 McCormick RoadCharlottesville, VA 22903
Office:B011+1-434-982-4729
Lab:B019+1-434-982-4751
Fax:+1-434-982-4766
WWW:http://www.people.virginia.edu/~mk9y/

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

Re: [R] Suppresing default text in pairs.lmList() in package = nlme

[Michael Kubovy]

Thanks Ken and Dieter,

I added xlab = '' and the text 'Scatter Plot Matrix' produced by pairs 
() applied to an lmList object (in nlme) went away. Skip the rest ---  
which is mainly autobiographical --- unless you're curious.

On Feb 11, 2007, at 4:25 AM, Dieter Menne wrote:

> Michael Kubovy  virginia.edu> writes:
>
>> I would like to suppress the text 'Scatter Plot Matrix' that appears
>> under the plot. Could someone please suggest how?
>
> Must be as special Virginia Brand of lmList. That Test does not  
> turn up in my
> output, and the only place I found it in the sources was as a comment
>
> ## scatter plot matrix plots, generally based on coef or ranef

Overnight I flew to New York (a slightly more left-leaning state),  
bought a new computer, installed R and all the packages, and then  
flew back with the new machine; the symptom remained  ;)

On Feb 11, 2007, at 5:20 AM, Ken Knoblauch wrote:

> It looks like you get this default xlab with splom so the title is
> misleading,
> though pairs.lmList calls splom in some cases, see near the end of the
> function
> where plotfun is assigned a value of either "xyplot" or "splom". I  
> was able
> to get the xlab to change by providing an explicit xlab = "machin"  
> argument
> to override the xlab argument in the call to splom and which is passed
> along by
> the ... argument.

Is it documented somewhere how to figure this out? It just didn't  
occur to me that the problem was as simple as changing the xlab and I  
don't yet have the habit of looking at code. Nevertheless I had done  
the following:

 > pairs
function (x, ...)
UseMethod("pairs")


... which did not make me happy. So then I tried:

 > pairs.lmList
Error: object "pairs.lmList" not found

... which did not make me happy. So then I tried:

 > apropos('pairs')
[1] "pairs" "pairs.default" "panel.pairs"
 > pairs.default
function (x, labels, panel = points, ..., lower.panel = panel,
...

... and then my eyes glazed over because nothing I saw seemed relevant.
_
Professor Michael Kubovy
University of Virginia
Department of Psychology
USPS: P.O.Box 400400Charlottesville, VA 22904-4400
Parcels:Room 102Gilmer Hall
 McCormick RoadCharlottesville, VA 22903
Office:B011+1-434-982-4729
Lab:B019+1-434-982-4751
Fax:+1-434-982-4766
WWW:http://www.people.virginia.edu/~mk9y/

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

[R] Suppresing default text in pairs.lmList() in package = nlme

[Michael Kubovy]

I would like to suppress the text 'Scatter Plot Matrix' that appears  
under the plot. Could someone please suggest how?
_
Professor Michael Kubovy
University of Virginia
Department of Psychology
USPS: P.O.Box 400400Charlottesville, VA 22904-4400
Parcels:Room 102Gilmer Hall
 McCormick RoadCharlottesville, VA 22903
Office:B011+1-434-982-4729
Lab:B019+1-434-982-4751
Fax:+1-434-982-4766
WWW:http://www.people.virginia.edu/~mk9y/

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

[R] OT: question about transformation

[Michael Kubovy]

If x is a r.v.
1/log(x) ~ N (approximately)
what does that tell me about the distribution of x (i.e., is it a  
standard distribution?)

Can you point me to a list where I would have appropriately posted this?
_
Professor Michael Kubovy
University of Virginia
Department of Psychology
USPS: P.O.Box 400400Charlottesville, VA 22904-4400
Parcels:Room 102Gilmer Hall
 McCormick RoadCharlottesville, VA 22903
Office:B011+1-434-982-4729
Lab:B019+1-434-982-4751
Fax:+1-434-982-4766
WWW:http://www.people.virginia.edu/~mk9y/

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

[R] How to add the variable name to a qqplot or densityplot in the diagonal of an splom?

[Michael Kubovy]

splom() doesn't complain here, but writes no names in the diagonal  
boxes. What am I missing?
I believe that I need to add something like grid.text(x, ...) to the  
diagonal panel, but I don't know how to get it cycle through the  
column labels. And should
 varname.col = 'blue', varname.cex = 1
be inside the diag.panel() function?

splom(szw[, n], pscales = 0,
 diag.panel = function(x, ...){
 panel.qqmathline(x, ...)
 panel.qqmath(x, ...)
 },
 lower.panel = function(x, y, ...){
 panel.xyplot(x, y, ..., col = 'lightblue')
 panel.loess(x, y, ..., col = 'red')
 },
 upper.panel = function(x, y, ...){
 panel.abline(lm(y~x),...)
 grid.text(round(cor(x, y, use = 'pairwise.complete.obs'), 2),
 x = unit(1, 'mm'),
 y = unit(1, 'npc') - unit(1, 'mm'),
 just = c('left', 'top'),
 gp = gpar(fontsize = 10))
     },
 varname.col = 'blue', varname.cex = 1
)


_
Professor Michael Kubovy
University of Virginia
Department of Psychology
USPS: P.O.Box 400400Charlottesville, VA 22904-4400
Parcels:Room 102Gilmer Hall
 McCormick RoadCharlottesville, VA 22903
Office:B011+1-434-982-4729
Lab:B019+1-434-982-4751
Fax:+1-434-982-4766
WWW:http://www.people.virginia.edu/~mk9y/

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

Re: [R] Replace individual values in a data frame with NA

[Michael Kubovy]

On Feb 9, 2007, at 11:10 AM, Jason Horn wrote:

> I'd like to replace a value in a data frame with an NA, but can't
> figure out how.
>
> For example, say you have
>
> a<-c(1,2,3,4)
> b<-c(5,6,7,8)
> data<-data.frame(a,b)
>
> Now, how would you set the third row of the second column ( data
> [[3,2]] ) to NA?

a<-c(1,2,3,4)
b<-c(5,6,7,8)
data<-data.frame(a,b)
data[3, 2] <- NA

?Extract
_
Professor Michael Kubovy
University of Virginia
Department of Psychology
USPS: P.O.Box 400400Charlottesville, VA 22904-4400
Parcels:Room 102Gilmer Hall
 McCormick RoadCharlottesville, VA 22903
Office:B011+1-434-982-4729
Lab:B019+1-434-982-4751
Fax:+1-434-982-4766
WWW:http://www.people.virginia.edu/~mk9y/

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

Re: [R] How to count the number of NAs in each column of a df?

[Michael Kubovy]

Dear Jim (25 minutes!), Richard (27 minutes!), and Chuck,

Thanks to your hints, I have come up with what I hope is a pithy  
idiom that drops columns of a dataframe (df) in which the number of  
NAs is > (e.g.) 30.

tmp <- df
tmp <- tmp[, which(as.numeric(colSums(is.na(tmp))) > 30)]
df <- tmp

I wonder if we have a place to keep R programming idioms (which  
probably get unnecessarily reinvented). Is the R-Wiki suitable?
_
Professor Michael Kubovy
University of Virginia
Department of Psychology
USPS: P.O.Box 400400Charlottesville, VA 22904-4400
Parcels:Room 102Gilmer Hall
 McCormick RoadCharlottesville, VA 22903
Office:B011+1-434-982-4729
Lab:B019+1-434-982-4751
Fax:+1-434-982-4766
WWW:http://www.people.virginia.edu/~mk9y/

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

[R] How to count the number of NAs in each column of a df?

[Michael Kubovy]

I would like to remove columns of a df which have too many NAs.

I think that summary() should give me the information, I just don't  
know how to access it.

Advice?
_
Professor Michael Kubovy
University of Virginia
Department of Psychology
USPS: P.O.Box 400400Charlottesville, VA 22904-4400
Parcels:Room 102Gilmer Hall
 McCormick RoadCharlottesville, VA 22903
Office:B011+1-434-982-4729
Lab:B019+1-434-982-4751
Fax:+1-434-982-4766
WWW:http://www.people.virginia.edu/~mk9y/

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

[R] Strange behavior of abc.ci() in boot

[Michael Kubovy]

require(boot)
# [1] TRUE
set.seed(1)
summary(rn <- rnorm(30))
# Min.  1st Qu.   Median Mean  3rd Qu. Max.
# -2.21500 -0.43500  0.25660  0.08246  0.70870  1.59500
abc.ci(rnorm(30), mean)
# [1] 0.9500 0.08245817 0.08245817

*
sessionInfo()
R version 2.4.1 (2006-12-18)
i386-apple-darwin8.8.1

locale:
C

attached base packages:
[1] "grid"  "datasets"  "stats" "graphics"  "grDevices"  
"utils" "methods"   "base"

other attached packages:
 boot   xtable latticeExtra  lattice  
gridBase MASS  JGR   iplots   JavaGD
 "1.2-27"  "1.4-3"  "0.1-4""0.14-16"  "0.4-3"  
"7.2-31" "1.4-15"  "1.0-5"  "0.3-5"
rJava
 "0.4-13"

_
Professor Michael Kubovy
University of Virginia
Department of Psychology
USPS: P.O.Box 400400Charlottesville, VA 22904-4400
Parcels:Room 102Gilmer Hall
 McCormick RoadCharlottesville, VA 22903
Office:B011+1-434-982-4729
Lab:B019+1-434-982-4751
Fax:+1-434-982-4766
WWW:http://www.people.virginia.edu/~mk9y/

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

[R] smartpred depends on fitted() in flexmix?

[Michael Kubovy]

Hi,

I was going through the examples in smartpred. It seems there's an  
unstated dependency on the fitted() function in package flexmix.

n = 20
set.seed(86)
x = sort(runif(n))
y = sort(runif(n))
library(splines)
fit = lm(y ~ ns(x, df=5))
plot(x, y)
lines(x, fitted(fit)) # won't work w/o prior loading of the flexmix  
package.
newx = seq(0, 1, len=n)
points(newx, predict(fit, data.frame(x=newx)), type="b", col=2, err=-1)
_____
Professor Michael Kubovy
University of Virginia
Department of Psychology
USPS: P.O.Box 400400Charlottesville, VA 22904-4400
Parcels:Room 102Gilmer Hall
 McCormick RoadCharlottesville, VA 22903
Office:B011+1-434-982-4729
Lab:B019+1-434-982-4751
Fax:+1-434-982-4766
WWW:http://www.people.virginia.edu/~mk9y/

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

Re: [R] lmer2 error under Mac OS X on PowerPC G5 but not on Dual-Core Intel Xeon

[Michael Kubovy]

On Jan 28, 2007, at 9:39 PM, Benilton Carvalho wrote:

> This seems to be due to the fact that you didn't have enough memory  
> when running lmer2.
>
> I might be wrong, but I think Calloc tries to get contiguous  
> memory, so this might the problem.
>
> If you are positive that you have enough memory, a gc() might help.

I have 2 GB memory on this machine. Should be enough, no?

 > gc()
   used (Mb) gc trigger (Mb) max used (Mb)
Ncells 1008175 27.01476915 39.5  1368491 36.6
Vcells  540055  4.21031040  7.9  1031026  7.9
 > (fm1 <- lmer2(Reaction ~ Days + (Days|Subject), sleepstudy))
Error in as.double(start) : Calloc could not allocate (903190944 of  
4) memory


> On Jan 28, 2007, at 8:35 PM, Michael Kubovy wrote:
>
>>> (fm1 <- lmer2(Reaction ~ Days + (Days|Subject), sleepstudy))
>> Error in as.double(start) : Calloc could not allocate (888475968 of
>> 4) memory
>> *
>>> sessionInfo()
>> R version 2.4.1 (2006-12-18)
>> powerpc-apple-darwin8.8.0
>>
>> locale:
>> C
>>
>> attached base packages:
>> [1] "grid"  "datasets"  "stats" "graphics"  "grDevices"
>> "utils" "methods"
>> [8] "base"
>>
>> other attached packages:
>>  lme4   Matrix   xtable latticeExtra  lattice
>> gridBase MASS
>> "0.9975-11"   "0.9975-8"  "1.4-3"  "0.1-4""0.14-16"
>> "0.4-3" "7.2-31"
>>   JGR   iplots   JavaGDrJava
>>  "1.4-15"  "1.0-5"  "0.3-5" "0.4-13"
>> *
>> lmer runs the example w/o a problem
>>
>> I just tried to run it on on Intel-based MacPro, and lmer2 ran
>> without a hitch.
_
Professor Michael Kubovy
University of Virginia
Department of Psychology
USPS: P.O.Box 400400Charlottesville, VA 22904-4400
Parcels:Room 102Gilmer Hall
 McCormick RoadCharlottesville, VA 22903
Office:B011+1-434-982-4729
Lab:B019+1-434-982-4751
Fax:+1-434-982-4766
WWW:http://www.people.virginia.edu/~mk9y/

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

[R] lmer2 error under Mac OS X on PowerPC G5 but not on Dual-Core Intel Xeon

[Michael Kubovy]

 > (fm1 <- lmer2(Reaction ~ Days + (Days|Subject), sleepstudy))
Error in as.double(start) : Calloc could not allocate (888475968 of  
4) memory
*
 > sessionInfo()
R version 2.4.1 (2006-12-18)
powerpc-apple-darwin8.8.0

locale:
C

attached base packages:
[1] "grid"  "datasets"  "stats" "graphics"  "grDevices"  
"utils" "methods"
[8] "base"

other attached packages:
 lme4   Matrix   xtable latticeExtra  lattice  
gridBase MASS
"0.9975-11"   "0.9975-8"  "1.4-3"  "0.1-4""0.14-16"   
"0.4-3" "7.2-31"
  JGR   iplots   JavaGDrJava
 "1.4-15"  "1.0-5"  "0.3-5" "0.4-13"
*
lmer runs the example w/o a problem

I just tried to run it on on Intel-based MacPro, and lmer2 ran  
without a hitch.
_
Professor Michael Kubovy
University of Virginia
Department of Psychology
USPS: P.O.Box 400400Charlottesville, VA 22904-4400
Parcels:Room 102Gilmer Hall
 McCormick RoadCharlottesville, VA 22903
Office:B011+1-434-982-4729
Lab:B019+1-434-982-4751
Fax:+1-434-982-4766
WWW:http://www.people.virginia.edu/~mk9y/

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

[R] Adding lines to xYplot

[Michael Kubovy]

I am using xYplot to plot data and CIs. How do I add several lines to  
the figure?
_
Professor Michael Kubovy
University of Virginia
Department of Psychology
USPS: P.O.Box 400400Charlottesville, VA 22904-4400
Parcels:Room 102Gilmer Hall
 McCormick RoadCharlottesville, VA 22903
Office:B011+1-434-982-4729
Lab:B019+1-434-982-4751
Fax:+1-434-982-4766
WWW:http://www.people.virginia.edu/~mk9y/

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

[R] Unbalanced design help

[Michael Kubovy]

Dear r-helpers,

I am looking at a designed experiment in which one predictor variable  
has 5 levels (0, ..., 4) and the other has 6 levels (1.1, ..., 1.6),  
with 33 observations per cell. This design was given to 13 subjects.
  0  1  2  3  4
   1.1 32 33  0  0  0
   1.2 33 33 33  0  0
   1.3 33 33 33 33  0
   1.4  0  0 33 33 33
   1.5  0  0  0 33 33
   1.6  0  0  0  0 33
The reason for this design is that low values of one predictor  
combined with high values of the other produce floor or ceiling  
effects at any reasonable sample size.

Can anyone point me to a theoretical discussion of the analysis of  
such data, and tools in R to use?
_
Professor Michael Kubovy
University of Virginia
Department of Psychology
USPS: P.O.Box 400400Charlottesville, VA 22904-4400
Parcels:Room 102Gilmer Hall
 McCormick RoadCharlottesville, VA 22903
Office:B011+1-434-982-4729
Lab:B019+1-434-982-4751
Fax:+1-434-982-4766
WWW:http://www.people.virginia.edu/~mk9y/



[[alternative HTML version deleted]]

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

Re: [R] Finding the effect of Box-Cox transformation using "vis.boxcoxu"

[Michael Kubovy]

?box.cox

?boxcox

On Jan 22, 2007, at 2:25 AM, Arun Kumar Saha wrote:

> I have a dataset 'data' and I want to see the effect of Box-Cox
> transformation on it Interactively for different lambda values. I  
> already
> got a look on function "vis.boxcoxu" in package "TeachingDemos". But I
> didn't find any option to put user's own dataset. Can anyone tell  
> me how to
> put my own dataset here i.e. "data"?

_
Professor Michael Kubovy
University of Virginia
Department of Psychology
USPS: P.O.Box 400400Charlottesville, VA 22904-4400
Parcels:Room 102Gilmer Hall
 McCormick RoadCharlottesville, VA 22903
Office:B011+1-434-982-4729
Lab:B019+1-434-982-4751
Fax:+1-434-982-4766
WWW:http://www.people.virginia.edu/~mk9y/

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

Re: [R] Newbie question: Statistical functions (e.g., mean, sd) in a "transform" statement?

[Michael Kubovy]

On Jan 19, 2007, at 12:54 PM, Ben Fairbank wrote:

> Given a data frame such as
>
> times
>
>time1time2 time3time4
> 1  70.408543 48.92378  7.399605 95.93050
> 2  17.231940 27.48530 82.962916 10.20619
> 3  20.279220 10.33575 66.209290 30.71846
> 4 NA 53.31993 12.398237 35.65782
> 5   9.295965   NA 48.929201   NA
> 6  63.966518 42.16304  1.777342   NA
>
> I cannot, however, find a way, other than "for" looping,
> to use statistical functions, such as mean or sd, to
> compute the new column.

times <- data.frame(time1 = rnorm(6, 50, 20), time2 = rnorm(6, 40, 15),
 time3 = rnorm(6, 60, 25), time4 = rnorm(6, 55, 23))
times[4,1] <- NA
times[5, c(2, 4)] <- NA
times[6, 4] <- NA
times$totaltime <- apply(times, 1, sum, na.rm = T)
times$meantime <- apply(times, 1, mean, na.rm = T)
times$sdtime <- apply(times, 1, sd, na.rm = T)

  time1time2time3time4 totaltime meantime   sdtime
1 28.84859 29.94037 92.11518 71.80472 222.70886 89.08354 71.11911
2 50.72260 39.02439 61.18364 31.63962 182.57024 73.02810 55.68944
3 11.75829 28.61262 72.37066 79.23817 191.97974 76.79189 62.99902
4   NA 27.23659 75.69952 38.19262 141.12872 70.56436 44.52787
5 31.05109   NA 52.41755   NA  83.46864 55.64576 21.52078
6 54.01291 52.48922 53.97689   NA 160.47902 80.23951 46.33038



_
Professor Michael Kubovy
University of Virginia
Department of Psychology
USPS: P.O.Box 400400Charlottesville, VA 22904-4400
Parcels:Room 102Gilmer Hall
 McCormick RoadCharlottesville, VA 22903
Office:B011+1-434-982-4729
Lab:B019+1-434-982-4751
Fax:+1-434-982-4766
WWW:http://www.people.virginia.edu/~mk9y/



[[alternative HTML version deleted]]

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

Re: [R] ANCOVA

[Michael Kubovy]

On Jan 6, 2007, at 8:34 AM, John Cardinale wrote:

> Are there any R function which can do analysis of covariance?

?lm
RSiteSearch('ancova')
_____
Professor Michael Kubovy
University of Virginia
Department of Psychology
USPS: P.O.Box 400400Charlottesville, VA 22904-4400
Parcels:Room 102Gilmer Hall
 McCormick RoadCharlottesville, VA 22903
Office:B011+1-434-982-4729
Lab:B019+1-434-982-4751
Fax:+1-434-982-4766
WWW:http://www.people.virginia.edu/~mk9y/

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

Re: [R] Help with filled.contour()

[Michael Kubovy]

To Jeff: a tip of the hat.

I have another question after Jeff's solution:

On Jan 3, 2007, at 5:09 PM, Jeffrey Horner wrote:

> Michael Kubovy  virginia.edu> writes:
>> I tried and it gave a strange result. See
>> http://people.virginia.edu/~mk9y/mySite/twoGaussian.R
>> and
>> http://people.virginia.edu/~mk9y/mySite/twoGaussian.pdf
>
>
> Try the image function. The smoothness of the plot will be  
> proportional to the length of x and y. For instance 200 isn't bad:
>
> mu1 <- 0
> mu2 <- 5
> s <- 1
> x <- seq(-2.5, 7.5, length = 200)
> y <- seq(-2.5, 2.5, length = 200)
> f <- function(x,y){
> term1 <- 1/(2*pi*sqrt(s*s))
> term2 <- -1/2
> term3 <- (x - mu1)^2/s
> term4 <- (y - mu1)^2/s
> term5 <- (x - mu2)^2/s
> term1*(.5 * exp(term2*(term3 + term4)) + .5 * exp(term2*(term5  
> + term4)))
> } # setting up the function of the multivariate normal density
> z <- outer(x, y, f)
> # persp(x, y, z)
> require(grDevices)
> #pdf('twoGaussian.pdf')
> #filled.contour(x, y, z, axes = F, frame.plot = F, asp = 1,
> #   col = gray(seq(0, 0.9, len = 25)), nlevels = 25)
> image(x,y,z,col=gray(seq(0,0.9,len=200)))

Is there a simpler way to get rid of axes, frame, and axis labels than
image(x, y, z, col = gray(seq(0, 0.9, len = 200)), asp = 1, xaxt =  
'n', yaxt = 'n', bty = 'n', ann = F)
?

_
Professor Michael Kubovy
University of Virginia
Department of Psychology
USPS: P.O.Box 400400Charlottesville, VA 22904-4400
Parcels:Room 102Gilmer Hall
 McCormick RoadCharlottesville, VA 22903
Office:B011+1-434-982-4729
Lab:B019+1-434-982-4751
Fax:+1-434-982-4766
WWW:http://www.people.virginia.edu/~mk9y/

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

Re: [R] slightly extended lm class

[Michael Kubovy]

On Jan 1, 2007, at 9:00 PM, ivo welch wrote:

> I have written a short lme.R function, which adds normalized
> coefficients and White heteroskedasticity-adjusted statistics to the
> standard output.  Otherwise, it behaves like lm.

Is it a good idea to call it lme, since there's a widely used lme()  
in the nlme package?
_____
Professor Michael Kubovy
University of Virginia
Department of Psychology
USPS: P.O.Box 400400Charlottesville, VA 22904-4400
Parcels:Room 102Gilmer Hall
 McCormick RoadCharlottesville, VA 22903
Office:B011+1-434-982-4729
Lab:B019+1-434-982-4751
Fax:+1-434-982-4766
WWW:http://www.people.virginia.edu/~mk9y/

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

Re: [R] Help with filled.contour()

[Michael Kubovy]

On Jan 1, 2007, at 2:13 PM, Dieter Menne wrote:

> Michael Kubovy  virginia.edu> writes:
>
>>
>> The following plot is a first approximation to what I need:
> .. see below
>
>> ***
>> It has four drawbacks. I would like:
>> (1) to remove the white grid;
>> (2) to remove the white contours;
>> (3) its appearance to be smooth;
>> (4) to remove the key.
>
>
> In partial fulfillment of your requirements, increasing the levels  
> gives a nice
> smooth plot.
>
> mu1 <- 0
> mu2 <- 5
> s <- 1
> x <- seq(-2.5, 7.5, length = 41)
> y <- seq(-2.5, 2.5, length = 41)
> f <- function(x,y){
> term1 <- 1/(2*pi*sqrt(s*s))
> term2 <- -1/2
> term3 <- (x - mu1)^2/s
> term4 <- (y - mu1)^2/s
> term5 <- (x - mu2)^2/s
> term1*(.5 * exp(term2*(term3 + term4)) + .5 * exp(term2*(term5 +
> term4)))
> }
>
> z <- outer(x, y, f)
> persp(x, y, z)
> require(grDevices)
> filled.contour(x, y, z, axes = F, frame.plot = F, asp = 1,
>  col = palette(gray(seq(0, 0.9, len = 500))), nlevels = 500)

Thanks, Dieter.

I tried and it gave a strange result. See
http://people.virginia.edu/~mk9y/mySite/twoGaussian.R
and
http://people.virginia.edu/~mk9y/mySite/twoGaussian.pdf

*
Session Info
*
 > sessionInfo()
R version 2.4.1 (2006-12-18)
powerpc-apple-darwin8.8.0

locale:
C

attached base packages:
[1] "datasets"  "grid"  "graphics"  "grDevices" "stats"  
"utils" "methods"
[8] "base"

other attached packages:
  JGR   iplots       JavaGD MASS  
gridBase  lattice
 "1.4-14"  "1.0-5"  "0.3-5" "7.2-30"  "0.4-3" 
"0.14-16"
latticeExtrarJava   xtable
  "0.1-4" "0.4-12"  "1.4-2"

_
Professor Michael Kubovy
University of Virginia
Department of Psychology
USPS: P.O.Box 400400Charlottesville, VA 22904-4400
Parcels:Room 102Gilmer Hall
 McCormick RoadCharlottesville, VA 22903
Office:B011+1-434-982-4729
Lab:B019+1-434-982-4751
Fax:+1-434-982-4766
WWW:http://www.people.virginia.edu/~mk9y/

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

[R] Help with filled.contour()

[Michael Kubovy]

The following plot is a first approximation to what I need:

***
mu1 <- 0
mu2 <- 5
s <- 1
x <- seq(-2.5, 7.5, length = 41)
y <- seq(-2.5, 2.5, length = 41)
f <- function(x,y){
term1 <- 1/(2*pi*sqrt(s*s))
term2 <- -1/2
term3 <- (x - mu1)^2/s
term4 <- (y - mu1)^2/s
term5 <- (x - mu2)^2/s
term1*(.5 * exp(term2*(term3 + term4)) + .5 * exp(term2*(term5 +  
term4)))
}
z <- outer(x, y, f)
persp(x, y, z)
require(grDevices)
filled.contour(x, y, z, axes = F, frame.plot = F, asp = 1,
 col = palette(gray(seq(0, 0.9, len = 25))), nlevels = 25)
***
It has four drawbacks. I would like:
(1) to remove the white grid;
(2) to remove the white contours;
(3) its appearance to be smooth;
(4) to remove the key.
_________
Professor Michael Kubovy
University of Virginia
Department of Psychology
USPS: P.O.Box 400400Charlottesville, VA 22904-4400
Parcels:Room 102Gilmer Hall
 McCormick RoadCharlottesville, VA 22903
Office:B011+1-434-982-4729
Lab:B019+1-434-982-4751
Fax:+1-434-982-4766
WWW:http://www.people.virginia.edu/~mk9y/

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

Re: [R] New to R

[Michael Kubovy]

On Dec 31, 2006, at 10:15 PM, Obinna Duru wrote:

> I have installed acepack but efforts to get started has been
> unsuccessful. I can't seem to be able to load my data files because I
> am yet to figure the syntax to use. Is there a work directory in R
> where I can put my files and call them anytime, like in Matlab? My
> files are on my C drive and I just can't figure the syntax to get  
> them into R.

Start by looking at
?getwd
?read.table
intall.packages('foreign') # if not installed
?read.spss
?read.ssd
etc.

and then
http://www.bioconductor.org/developers/progRef/fileHandling.pdf

followed by the more general
http://cran.r-project.org/doc/manuals/R-intro.pdf
_________
Professor Michael Kubovy
University of Virginia
Department of Psychology
USPS: P.O.Box 400400Charlottesville, VA 22904-4400
Parcels:Room 102Gilmer Hall
 McCormick RoadCharlottesville, VA 22903
Office:B011+1-434-982-4729
Lab:B019+1-434-982-4751
Fax:+1-434-982-4766
WWW:http://www.people.virginia.edu/~mk9y/

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

Re: [R] Difference of array and vector

[Michael Kubovy]

On Dec 29, 2006, at 2:16 PM, Geoffrey Zhu wrote:

> Does anyone know what is the real difference of vector and array?

# Here's a vector:
letv <- letters
is.vector(letv)

# Here's an array:
leta <- as.array(letters)
is.array(leta)
attributes(letv)
attributes(leta)

To understand the importance of attributes:
http://www.burns-stat.com/pages/Spoetry/essentials.pdf
There Burns writes (p. 3):
"An S array is merely a vector that has a dim attribute, and  
optionally a dimnames attribute."
_________
Professor Michael Kubovy
University of Virginia
Department of Psychology
USPS: P.O.Box 400400Charlottesville, VA 22904-4400
Parcels:Room 102Gilmer Hall
 McCormick RoadCharlottesville, VA 22903
Office:B011+1-434-982-4729
Lab:B019+1-434-982-4751
Fax:+1-434-982-4766
WWW:http://www.people.virginia.edu/~mk9y/



[[alternative HTML version deleted]]

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

Re: [R] Rating competitors

[Michael Kubovy]

I would start with elimination-by-aspects models:
?eba
I would read the Tversky 1972 paper (cited on the help page for the  
eba() function), which is brilliant.

Jeff Newmiller wrote:

> I am looking for hints on how to estimate ratings for competitors
> in an ongoing pairwise competition using R... my particular area of
> interest being the game of Go, but the idea of identifying ratings
> (on a continuous scale) rather than relative rankings seems easily
> generalized to other competitions so I thought someone might be
> studying something related already.
>
> I presume the rating of a competitor would be best modeled as a random
> variate on the rating scale, and an encounter between two
> competitors would be represented by a binary result.  Logistic  
> regression
> seems promising, but I am at a loss how to represent the model since
> the pairings are arbitrary and not necessarily repeated often.
_________
Professor Michael Kubovy
University of Virginia
Department of Psychology
USPS: P.O.Box 400400Charlottesville, VA 22904-4400
Parcels:Room 102Gilmer Hall
 McCormick RoadCharlottesville, VA 22903
Office:B011+1-434-982-4729
Lab:B019+1-434-982-4751
Fax:+1-434-982-4766
WWW:http://www.people.virginia.edu/~mk9y/

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

Re: [R] sequential row selection in dataframe

[Michael Kubovy]

On Dec 26, 2006, at 12:07 AM, Pedro Mardones wrote:

> I'm wondering if there is any 'efficient' approach for selecting a
> sample of 'every nth rows'  from a dataframe. For example, let's use
> the dataframe GAGurine in MASS library:
>
>> length(GAGurine[,1])
> [1] 314
>
> # select an 75% of the dataset, i.e. = 236 rows, every 2 rows starting
> from row 1
>> test<-GAGurine[seq(1,314,2),]
>> length(test[,1])
> [1] 157
>
> # so, I still need another 79 rows, one way could be:
> test2<-GAGurine[-seq(1,314,2),]
>> length(test2[,1])
> [1] 157
>> test3<-test2[seq(1,157,2),]
>
> # and then
> final<-rbind(test2,test3)
>> length(final[,1])
> [1] 236
>
> Does anyone have a better idea to get the same results but without
> creating different datasets like test2 and test3?

A probabilistic approach:

len <- length(GAGurine[,1])
GAGu <- GAGurine[sample(1:len, round(.75 * len)), ] # 236 rows

A deterministic one:

nr <- 1 #or 2
GAGu2 <- GAGurine[-seq(nr, len, 4),] # drop every 4th, giving 235 rows
nr <- 3 # or 4
will give 236 rows.
_
Professor Michael Kubovy
University of Virginia
Department of Psychology
USPS: P.O.Box 400400Charlottesville, VA 22904-4400
Parcels:Room 102Gilmer Hall
 McCormick RoadCharlottesville, VA 22903
Office:B011+1-434-982-4729
Lab:B019+1-434-982-4751
Fax:+1-434-982-4766
WWW:http://www.people.virginia.edu/~mk9y/

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

Re: [R] complex barplot enquiry

[Michael Kubovy]

On Dec 23, 2006, at 12:07 AM, Bob Green wrote:



> The data consists of a series of pre and post variables, in a  
> dataframe
> called 'offend'. I am interested in graphically depicting the pre &  
> post
> values for a factor variable called 'decision' which has 4 values :  
> nusm,
> fit, unsound & unfit. An example of a pre and post variable is:  
> pre.damage
> and post.damage (data below). I wanted to exclude all values = 0, and
> collapse the remaining values into 3 categories .e.g 1, 2-5, 6 plus.
>
> For each 'decision' group the final graph would contain the  
> alternating pre
> and post values for 'damage'.  For example, for 'fit' the first bar  
> would
> have a pre value = 18 & a post value = 4, the second pre bar would  
> = 17 and
> the post bar = 2, and the third pre bar would = 2  & the post = 1.



> > table(decision,pre.damage)
>   pre.damage
> decision0  1   2   3   4   5   6   7   8   9  10  11  13   
> 14 20
>fit   89   189   4   2   2   1   0   0   1   0   0
> 0   0   0
>nusm 158  44  15   5   5   6   6   3   1   1   2   0   0
> 0   0
>unfit 272   2   0   0   1   0   0   0   0   0   0
> 0   0   0
>unsound 333  68  36   8   5   4   1   0   0   0   0   1   1   1   1
>
> > table(decision,post.damage)
>   post.damage
> decision0  12   3   4   5   6   7  10  12  13
>fit 11941   1   0   0   0   0   1   0   0
>nusm 212  16   11   2   0   2   1   1   0   1   0
>unfit 3110   0   0   0   0   0   0   0   0
>unsound 44185   3   1   0   0   0   0   0   1

Could you please post commands to define the df 'offend'?
_
Professor Michael Kubovy
University of Virginia
Department of Psychology
USPS: P.O.Box 400400Charlottesville, VA 22904-4400
Parcels:Room 102Gilmer Hall
 McCormick RoadCharlottesville, VA 22903
Office:B011+1-434-982-4729
Lab:B019+1-434-982-4751
Fax:+1-434-982-4766
WWW:http://www.people.virginia.edu/~mk9y/

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

Re: [R] Newbie data organisation/structures question...

[Michael Kubovy]

On Dec 20, 2006, at 11:05 AM, Gav Wood wrote:

> So my data is in this sort of format:
>
> P  T  I
> 1  1  (1, 2, 3)
> 2  1  (2, 4)
> 1  2  (1, 3, 6, 7)
> 2  2  (6)

Not knowing why you organized the data as you did, let me suggest  
another approach:

iv <- c(1, 2, 3, 2, 4, 1, 3, 6, 7, 6)
p <- c(1, 1, 1, 2, 2, 1, 1, 1, 1, 2)
t <- rep(1:2, each = 5)
dat <- data.frame(iv, p, t)

> And I want to be able to quickly get:
>
> The I when both P and T are given. e.g.:
> P = 2, T = 2; I = (6)

subset(dat, p == 2 & t ==2)$iv

> The concatenated vector of Is when P and a subset of T is given, e.g.:
> P = 1, T = 1:2;  Is = (1, 2, 3, 1, 3, 6, 7)

(iv1 <- subset(dat, p == 1)$iv)

> The length of that vector.

length(iv1)

> A list of Is when either P or T is given. e.g.:
> P = 2: I = (2, 4), (6)
> T = 1: I = (1, 2, 3), (1, 3, 6, 7)

list(p2t1 = subset(dat, p == 2 & t ==1)$iv, p2t2 = subset(dat, p == 2  
& t ==2)$iv)
list(p1t1 = subset(dat, p == 1 & t ==1)$iv, p1t2 = subset(dat, p == 1  
& t ==2)$iv) # correcting your requirement to get your result

There are many other ways of getting the results you need as Marc  
Schwartz pointed out in his reply.
_
Professor Michael Kubovy
University of Virginia
Department of Psychology
USPS: P.O.Box 400400Charlottesville, VA 22904-4400
Parcels:Room 102Gilmer Hall
 McCormick RoadCharlottesville, VA 22903
Office:B011+1-434-982-4729
Lab:B019+1-434-982-4751
Fax:+1-434-982-4766
WWW:http://www.people.virginia.edu/~mk9y/

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

[R] Problem with glmmADMB

[Michael Kubovy]

library(glmmADMB)
#Example for   glmm.admb
data(epil2)
glmm.admb(y~Base*trt+Age 
+Visit,random=~Visit,group="subject",data=epil2,family="nbinom")

Gives:
Error in glmm.admb(y ~ Base * trt + Age + Visit, random = ~Visit,  
group = "subject",  :
The function maximizer failed

**
R version 2.4.1 RC (2006-12-14 r40181)
powerpc-apple-darwin8.8.0

locale:
C

attached base packages:
[1] "utils" "stats" "graphics"  "grDevices" "datasets"   
"methods"   "base"

other attached packages:
glmmADMB   JGRiplotsJavaGD rJava  MASS   lattice
 "0.3"  "1.4-14"   "1.0-5"   "0.3-5"  "0.4-12"  "7.2-30" "0.14-16"

_
Professor Michael Kubovy
University of Virginia
Department of Psychology
USPS: P.O.Box 400400Charlottesville, VA 22904-4400
Parcels:Room 102Gilmer Hall
 McCormick RoadCharlottesville, VA 22903
Office:B011+1-434-982-4729
Lab:B019+1-434-982-4751
Fax:+1-434-982-4766
WWW:http://www.people.virginia.edu/~mk9y/

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

Re: [R] how to replace some objects?

[Michael Kubovy]

On Dec 19, 2006, at 3:05 AM, Zhang Jian wrote:

> I want to replace some objects in one row or column.For example,
> One colume: a,b,a,c,b,b,a,a,c.
> I want to replace "a" with "1", "b" with "2", and "c" with "3".
> Like this: 1,2,1,3,2,2,1,1,3.

let <- c('a', 'b', 'a', 'c', 'b', 'b', 'a', 'a', 'c')
library(car)
num <- recode(let, " 'a' = 1; 'b' = 2; else = 3 ")

_
Professor Michael Kubovy
University of Virginia
Department of Psychology
USPS: P.O.Box 400400Charlottesville, VA 22904-4400
Parcels:Room 102Gilmer Hall
 McCormick RoadCharlottesville, VA 22903
Office:B011+1-434-982-4729
Lab:B019+1-434-982-4751
Fax:+1-434-982-4766
WWW:http://www.people.virginia.edu/~mk9y/

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

Re: [R] digital length

[Michael Kubovy]

?Round

On Dec 14, 2006, at 2:52 AM, XinMeng wrote:

> How can I control the digital length of data?
>
> e.g:
> 0.1234 is the output of an algorithm.
> What I want is 0.12 instead.

_
Professor Michael Kubovy
University of Virginia
Department of Psychology
USPS: P.O.Box 400400Charlottesville, VA 22904-4400
Parcels:Room 102Gilmer Hall
 McCormick RoadCharlottesville, VA 22903
Office:B011+1-434-982-4729
Lab:B019+1-434-982-4751
Fax:+1-434-982-4766
WWW:http://www.people.virginia.edu/~mk9y/

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

[R] Sweave and warning messages

[Michael Kubovy]

How does one tell Sweave() to include analysis warning messages in  
the verbatim output?
_
Professor Michael Kubovy
University of Virginia
Department of Psychology
USPS: P.O.Box 400400Charlottesville, VA 22904-4400
Parcels:Room 102Gilmer Hall
 McCormick RoadCharlottesville, VA 22903
Office:B011+1-434-982-4729
Lab:B019+1-434-982-4751
Fax:+1-434-982-4766
WWW:http://www.people.virginia.edu/~mk9y/

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

Re: [R] A smal fitting problem...

[Michael Kubovy]

On Dec 8, 2006, at 7:42 AM, David Barron wrote:

> b = mean(y)
>
> On 08/12/06, Kåre Edvardsen <[EMAIL PROTECTED]> wrote:
>> Dear R-helpers,
>>
>> I'm for sure not familiar with R, but it seem like a nice sofware  
>> tool,
>> so I've decided to try using it.
>>
>> Here is my problem I just can't figure out:
>>
>> I'd like to do least square fit of a straight horizontal (a = 0)  
>> line y
>> = ax + b through some data points
>>
>> x = (3,4,5,6,7,8)
>>
>> y = (0.62, 0.99, 0.83, 0.69, 0.76, 0.82)
>>
>> How would i find b?

And in the context of linear models:
x <- 3:8
y <- c(0.62, 0.99, 0.83, 0.69, 0.76, 0.82)
lm(y ~ 1)
_
Professor Michael Kubovy
University of Virginia
Department of Psychology
USPS: P.O.Box 400400Charlottesville, VA 22904-4400
Parcels:Room 102Gilmer Hall
 McCormick RoadCharlottesville, VA 22903
Office:B011+1-434-982-4729
Lab:B019+1-434-982-4751
Fax:+1-434-982-4766
WWW:http://www.people.virginia.edu/~mk9y/

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

Re: [R] Help with response CIs for lme

[Michael Kubovy]

Thanks, Spencer, for spending time on my question.

On Dec 5, 2006, at 1:48 AM, Spencer Graves wrote:

>  No, your example computation does NOT produce the desired "lme  
> prediction intervals".  I ran your script and got exactly the same  
> numbers for upper and lower limits.  Even without any consideration  
> of statistical theory, this suggests either shockingly precise  
> statistical estimation or a problem with your algorithm.

There may have been a typo in what I sent to the list.
The following script produces prediction intervals that are not  
nearly as shocking as intervals with 0 width:

library(mlmRev)
library(nlme)
hsb.lme <- lme(mAch ~ minrty * sector, random = ~ 1 | cses, Hsb82)
(hsb.int <- intervals(hsb.lme))
hsb.new <- data.frame(
  minrty = rep(c('No', 'Yes'), 2),
  sector = rep(c('Public', 'Catholic'), each = 2))
cbind(hsb.new, predict(hsb.lme, hsb.new, level = 0))

cbind(hsb.new, rbind(
 hsb.int[[1]][1,],
 hsb.int[[1]][1,] + hsb.int[[1]][2,],
 hsb.int[[1]][1,] + hsb.int[[1]][3,],
 hsb.int[[1]][1,] + hsb.int[[1]][2,] + hsb.int[[1]][3,] + hsb.int 
[[1]][4,]
))

The latter produces:
   minrty   sector lower  est.upper
1 No   Public 12.467015 12.712155 12.95730
2Yes   Public  6.657422  7.387386  8.11735
3 No Catholic 14.54 15.123113 15.72400
4Yes Catholic 10.062856 11.821924 13.58099



>  To briefly outline some of the difficulties, we first should  
> ask if you want confidence intervals for the MEAN of future values  
> or for the future values themselves?

I want the former, and that's what I think I computed.

> And how do you want to handle the random effects?  If you want the  
> mean of the fixed effects, averaging over random effects, that  
> should be fairly easy.

Exactly. And my question is --- if the CIs I obtained were correct  
--- is there a less laborious way to obtain them?
_
Professor Michael Kubovy
University of Virginia
Department of Psychology
USPS: P.O.Box 400400Charlottesville, VA 22904-4400
Parcels:Room 102Gilmer Hall
 McCormick RoadCharlottesville, VA 22903
Office:B011+1-434-982-4729
Lab:B019+1-434-982-4751
Fax:+1-434-982-4766
WWW:http://www.people.virginia.edu/~mk9y/



[[alternative HTML version deleted]]

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

Re: [R] stat question - not R question so ignore if not interested

[Michael Kubovy]

On Dec 5, 2006, at 3:42 PM, Leeds, Mark ((IED)) wrote:

> If do a scattrplot of data ( x and y ) and there are two clouds of
> points. One cloud is in the left
> bottom corner of the plot and the other cloud is in the upper right.
>
> If I fit a regression line to this data ( or equivalently ,  
> calculate a
> correlation ), then obviously, it is going to seem like
> x and y are related because a line has to be connected between the 2
> clouds. But, there must be a regression assumption that
> is violated here because if the regressions are done separately on  
> each
> cloud, then there really isn't
> a relationship between x and y. I was just wondering 1) what  
> assumption
> in regression is being violated in
> the first case or 2) possibly if the regression is valid and the  
> results
> just have some different interpreation ?

One needs only to look at diagnostic plots:

Suppose
set.seed(2)
xy <- data.frame(y = c(rnorm(300), rnorm(300, 5)), x = c(rnorm(300),  
rnorm(300, 5)))
op <- par(mfrow = c(2,2))
plot(lm(y ~ x, xy))
par(op)

The model does not fit well because the residuals aren't flat as a  
function of fit and because homoscedasticity is violated.

When this happens we might try a different approach:
require(sm)
xy.sm <- sm.regression(xy$x, xy$y)

Whenever there's a big discrepancy between an OLS fit and a robust  
one, we should not pursue the OLS one w/o reinterpretation, which  
others have discussed in their replies.
_
Professor Michael Kubovy
University of Virginia
Department of Psychology
USPS: P.O.Box 400400Charlottesville, VA 22904-4400
Parcels:Room 102Gilmer Hall
 McCormick RoadCharlottesville, VA 22903
Office:B011+1-434-982-4729
Lab:B019+1-434-982-4751
Fax:+1-434-982-4766
WWW:http://www.people.virginia.edu/~mk9y/

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

[R] overplot() examples fail (package:gplots)

[Michael Kubovy]

require(gplots)
data(rtPCR)

overplot( RQ ~ Conc..ug.ml. | Test.Substance,
data=rtPCR,
subset=Detector=="ProbeType 7" & Conc..ug.ml. > 0,
same.scale=TRUE, log="xy", f=3/4,
main="Detector=ProbeType 7",
xlab="Concentration (ug/ml)",
ylab="Relative Gene Quantification"
)

# Error in lowess.default(mf[[-response]], mf[[response]], f = f,  
iter = iter,  :
#   'x' is empty
# In addition: Warning messages:
# 1: no non-missing arguments to min; returning Inf
# 2: no non-missing arguments to max; returning -Inf

**
# The second example fails as well, but with:
# Error in plot.window(xlim, ylim, log, asp, ...) :
#   need finite 'ylim' values
# In addition: Warning messages:
# 1: no non-missing arguments to min; returning Inf
# 2: no non-missing arguments to max; returning -Inf

**
R version 2.4.0 (2006-10-03)
powerpc-apple-darwin8.7.0

locale:
C

attached base packages:
[1] "splines"   "grid"  "methods"   "utils" "stats"  
"graphics"  "grDevices" "datasets"
[9] "base"

other attached packages:
 gplots  gdata gtools   lme4gam 
effects  Hmisc   codagmodels
"2.3.2""2.3.1""2.3.0" "0.9975-9" "0.97""1.0-9" 
"3.1-2"   "0.10-7"   "2.13.1"
 MatrixreshapeJGR iplots JavaGD   
rJava   MASSlattice
"0.9975-6""0.7.1"   "1.4-13""1.0-5""0.3-5"   "0.4-12"
"7.2-30"  "0.14-15"

_
Professor Michael Kubovy
University of Virginia
Department of Psychology
USPS: P.O.Box 400400Charlottesville, VA 22904-4400
Parcels:Room 102Gilmer Hall
 McCormick RoadCharlottesville, VA 22903
Office:B011+1-434-982-4729
Lab:B019+1-434-982-4751
Fax:+1-434-982-4766
WWW:http://www.people.virginia.edu/~mk9y/

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

[R] Problem with pairs() in nlme

[Michael Kubovy]

Dear r-helpers,

After successfully running
require(nlme)
vfr.lmL <- lmList(
 estimate ~ (slant + respType + visField + hand)^2 | subject, vfr
 )

pairs(vfr.lmL, id = 0.01, adj = -0.5) # Pinheiro & Bates (p. 141)  
produces the following error:
Error in sprintf(gettext(fmt, domain = domain), ...) :
object "form" not found

Any guesses as to what I may have done wrong?

vfr is not a grouped object because:

vfrg <- groupedData(estimate ~ slant | subject, data = vfr)
vfrg.lmL <- lmList(
estimate ~ (slant + respType + visField + hand)^2 | subject,
vfrg)
Error in function (classes, fdef, mtable)  :
unable to find an inherited method for function "lmList", for  
signature "formula", "nfnGroupedData"

The data are at
http://www.people.virginia.edu/~mk9y/LongFormat.csv

To get the working data.frame:
vf <- read.table("http://www.people.virginia.edu/~mk9y/ 
LongFormat.csv", header = TRUE, sep=',')
names(vf) <- c(
 'slant', 'estimate', 'respType', 'visField', 'hand',
 'order', 'sex', 'subject', 'handPref'
)
require(reshape)
vf <- sort_df(vf, vars = 'slant')
vf$subject <- factor(vf$subject)
summary(vf)
vfr <- subset(vf, handPref == 'R', select = -handPref)
vfr <- subset(vfr, slant %in% c(10, 20, 30))

_
Professor Michael Kubovy
University of Virginia
Department of Psychology
USPS: P.O.Box 400400Charlottesville, VA 22904-4400
Parcels:Room 102Gilmer Hall
 McCormick RoadCharlottesville, VA 22903
Office:B011+1-434-982-4729
Lab:B019+1-434-982-4751
Fax:+1-434-982-4766
WWW:http://www.people.virginia.edu/~mk9y/

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

[R] Slight discrepancy between predict.lm() and all.effects()

[Michael Kubovy]

In the course of exploring response prediction, I stumbled upon a  
small discrepancy between the CIs produced by predict.lm() and  
all.effects()

require(mlmRev)
require(effects)
hsb.lm <- lm(mAch ~ minrty * sector, Hsb82)
hsb.new <- data.frame(
 minrty = rep(c('No', 'Yes'), 2),
 sector = rep(c('Public', 'Catholic'), each = 2))
hsb.eff <- all.effects(hsb.lm)
cbind(
 hsb.new,
 predict(hsb.lm, hsb.new, interval = 'confidence', type =  
'response')
)
# the following lower and upper bounds differ starting with the  
fourth decimal place
data.frame(
 hsb.new,
 fit = hsb.eff[[1]]$fit,
 lwr = hsb.eff[[1]]$lower,
 upr = hsb.eff[[1]]$upper
)

Is this due to rounding or algorithm?
_
Professor Michael Kubovy
University of Virginia
Department of Psychology
USPS: P.O.Box 400400Charlottesville, VA 22904-4400
Parcels:Room 102Gilmer Hall
 McCormick RoadCharlottesville, VA 22903
Office:B011+1-434-982-4729
Lab:B019+1-434-982-4751
Fax:+1-434-982-4766
WWW:http://www.people.virginia.edu/~mk9y/

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

[R] Help with response CIs for lme

[Michael Kubovy]

Hi,

Can someone please offer a procedure for going from CIs produced by  
intervals.lme() to fixed-effects response CIs.

Here's a simple example:

library(mlmRev)
library(nlme)
hsb.lme <- lme(mAch ~ minrty * sector, random = ~ 1 | cses, Hsb82)
(intervals(hsb.lme))
(hsb.new <- data.frame
 minrty = rep(c('No', 'Yes'), 2),
 sector = rep(c('Public', 'Catholic'), each = 2)))
cbind(hsb.new, predict(hsb.lme, hsb.new, level = 0))

Is the following correct (I know from the previous command that the  
estimate is correct)?
cbind(hsb.new, rbind(hsb.int[[1]][1,], hsb.int[[1]][1,]+hsb.int[[1]] 
[2,], hsb.int[[1]][1,]+hsb.int[[1]][3,], hsb.int[[1]][1,]+hsb.int[[1]] 
[2,] + hsb.int[[1]][3,] + hsb.int[[1]][4,]))

If so, is there an easier way to write it?
_____
Professor Michael Kubovy
University of Virginia
Department of Psychology
USPS: P.O.Box 400400Charlottesville, VA 22904-4400
Parcels:Room 102Gilmer Hall
 McCormick RoadCharlottesville, VA 22903
Office:B011+1-434-982-4729
Lab:B019+1-434-982-4751
Fax:+1-434-982-4766
WWW:http://www.people.virginia.edu/~mk9y/

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

[R] Can plot(augPred()) be used with objects that are not groupedData?

[Michael Kubovy]

Dear Friends,

Is there a workaround to this problem or must the object plotted by  
plot(augPred(object, ...)) be a grouped object?

score <- c(108, 103, 96, 84, 118, 110, 129, 90, 84, 96,
 105, 113, 96, 117, 107, 85, 125, 107, 128, 84, 104,
 100, 114, 117, 110, 127, 106, 92, 125, 96, 123, 101,
 100, 103, 105, 132, 122, 133, 107, 99, 116, 91, 128,
 113, 88, 105, 112, 130)
mcc <- data.frame(
subject = factor(rep(1:12, 4)),
age =rep(seq(30, 48, by = 6), each = 12),
score = score
)
mccG <- groupedData(score ~ age | subject, data = mcc)
class(mcc)
class(mccG)
mcc.lme <- lme(score ~ age, random = ~ 1 | subject, mcc)
mccG.lme <- lme(mccG)
summary(augPred(mcc.lme, primary = ~ mcc$age)) # works
summary(augPred(mccG.lme)) # works
plot(augPred(mcc.lme, primary = ~ mcc$age)) # doesn't show fitted line
plot(augPred(mccG.lme)) # does show fitted line

_____
Professor Michael Kubovy
University of Virginia
Department of Psychology
USPS: P.O.Box 400400Charlottesville, VA 22904-4400
Parcels:Room 102Gilmer Hall
 McCormick RoadCharlottesville, VA 22903
Office:B011+1-434-982-4729
Lab:B019+1-434-982-4751
Fax:+1-434-982-4766
WWW:http://www.people.virginia.edu/~mk9y/

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

[R] Help with plot(augPred()) plot(comparePred()) in nlme

[Michael Kubovy]

Dear r-helpers,

When I issue the commands
plot(augPred(mcc.lme, primary = ~ mcc$age, length.out = 2, level = c 
(0, 1)))
or
plot(comparePred(mcc1.lis, mcc.lme, primary = ~ mcc$age, length.out =  
2), layout = c(4, 3)))
no lines are drawn.

But
fm1 <- lme(Orthodont)
plot(augPred(fm1, level = 0:1, length.out = 2))
does plot lines.

augPred() and comparePredict() seem to be working correctly, as shown  
below

augPred(mcc.lme, primary = ~ mcc$age, length.out = 2, level = c(0, 1))
   mcc$age .groups score.type
1  30   1   108 original
2  30   2   103 original
3  30   396 original
4  30   484 original
5  30   5   118 original
6  30   6   110 original
...
91  30  10   101 predict.subject
92  48  10   101 predict.subject
93  30  11   109 predict.subject
94  48  11   109 predict.subject
95  30  12   122 predict.subject
96  48  12   122 predict.subject

comparePred(mcc1.lis, mcc.lme, primary = ~ mcc$age, length.out = 2))
   mcc$age subject score.type
1  30   1   108 original
2  30   296 original
3  30   3   110 original
4  30   4   122 original
5  30   5   103 original
6  30   6   117 original
...
671  30  10   101 mcc.lme
681  48  10   101 mcc.lme
691  30  11   109 mcc.lme
701  48  11   109 mcc.lme
711  30  12   122 mcc.lme
721  48  12   122 mcc.lme


_
Professor Michael Kubovy
University of Virginia
Department of Psychology
USPS: P.O.Box 400400Charlottesville, VA 22904-4400
Parcels:Room 102Gilmer Hall
 McCormick RoadCharlottesville, VA 22903
Office:B011+1-434-982-4729
Lab:B019+1-434-982-4751
Fax:+1-434-982-4766
WWW:http://www.people.virginia.edu/~mk9y/

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

Re: [R] have I an actual matrix?

[Michael Kubovy]

On Nov 7, 2006, at 6:25 PM, Ricardo Rodríguez - Your EPEC ICT Team  
wrote:

> library(RMySQL)
> con <- dbConnect(dbDriver
> ("MySQL"),host='localhost',username='root',dbname='ibdona')
> rs <- dbGetQuery (con,"select n,year from ibdona.library_location")
> dbDisconnect(con)
> Graph <- barplot(rs)
>
> And here the error I get...
>
> Error in barplot.default(rs) : 'height' must be a vector or a matrix
>
> paste(rs) gives me...
>
> [1] "c(307, 65, 2, 28, 3, 229, 81, 5, 7, 558, 134, 53, 9)"
> [2] "c(2002, 2002, 2002, 2002, 2002, 2003, 2003, 2003, 2003, 2004,
> 2004, 2004, 2004)"

nums <- c(307, 65, 2, 28, 3, 229, 81, 5, 7, 558, 134, 53, 9)
names(nums) <- c(2002, 2002, 2002, 2002, 2002, 2003, 2003, 2003,  
2003, 2004, 2004, 2004, 2004)
barplot(nums)

or

require(gplots)
barplot2(nums, plot.grid = TRUE, las = 1) # offers many nice options
_
Professor Michael Kubovy
University of Virginia
Department of Psychology
USPS: P.O.Box 400400Charlottesville, VA 22904-4400
Parcels:Room 102Gilmer Hall
 McCormick RoadCharlottesville, VA 22903
Office:B011+1-434-982-4729
Lab:B019+1-434-982-4751
Fax:+1-434-982-4766
WWW:http://www.people.virginia.edu/~mk9y/

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

Re: [R] Equivalent to Stata command

[Michael Kubovy]

On Nov 4, 2006, at 10:57 AM, Jon Minton wrote:

> Do any of you know if there's an equivalent function to the Stata  
> command:
>
> Use v1 v3 v5 v2 v7 using file.dat
>
> In R? (i.e. something that just selectively loads certain variables  
> from a
> file and not others)?

How about (assuming that columns in 'file.dat' are separated by white  
space)
data <- read.table("file.dat", header = T, sep = " ") # reads in the  
whole file
data <- data[, c(1, 3, 5, 2, 7)] # keeps the columns you need

For more, see http://finzi.psych.upenn.edu/R/doc/manual/R-data.html
_
Professor Michael Kubovy
University of Virginia
Department of Psychology
USPS: P.O.Box 400400Charlottesville, VA 22904-4400
Parcels:Room 102Gilmer Hall
 McCormick RoadCharlottesville, VA 22903
Office:B011+1-434-982-4729
Lab:B019+1-434-982-4751
Fax:+1-434-982-4766
WWW:http://www.people.virginia.edu/~mk9y/

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

[R] CI -- fixed effects -- response -- mixed models

[Michael Kubovy]

Dear Friends,

I have been following the discussions of mixed-effects models with  
some confusion, and I realize that much of this is work in progress,  
and that much of the discussion is beyond my knowledge of statistics.

My question is simple, though: Is there a set of commands that will  
produce an output equivalent to the exceedingly useful
predict(bl.lm, newdata = bl.new, type = 'response', interval =  
'confidence')
for (at least) gaussian mixed-effects model (but preferably for glmm  
models)? I have looked at
predict.glmmPQL {MASS}
predict.lme {nlme}
and neither of them offer exactly what I need, and I just don't know  
enough to go from what they do offer to what I need. I have Pinheiro  
& Bates as well as MASS, and both come tantalizingly close to what I  
need, but I can't figure out the next step.

AFAIK, there's nothing of the sort for aov(. ~ ., Error(...) ...)  
either.

In order to publish results of designed psychological experiments  
(most of which are of the classic anova variety, with all the  
predictors being factors), we need to plot error bars on our  
interaction plots. I suspect (from the discussion on this list) that  
my colleagues are using SPSS or SAS, and reporting incorrect CIs.

Not being a statistician I have hit a wall here. I'm not sure if the  
transition to the results I need is staring me in the face, and I  
don't know enough to take the next step, or if the tools aren't yet  
available. In any event, I would be very grateful to you for guidance  
on how to proceed.
_
Professor Michael Kubovy
University of Virginia
Department of Psychology
USPS: P.O.Box 400400Charlottesville, VA 22904-4400
Parcels:Room 102Gilmer Hall
 McCormick RoadCharlottesville, VA 22903
Office:B011+1-434-982-4729
Lab:B019+1-434-982-4751
Fax:+1-434-982-4766
WWW:http://www.people.virginia.edu/~mk9y/

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

Re: [R] no possible to load a package correctly dowloaded by the R prompt

[Michael Kubovy]

Hi,

On Oct 30, 2006, at 11:43 AM, Tancredi Caruso wrote:

>> install.packages(c(untb))
>
>> library(untb)

> Carico il pacchetto richiesto: "partitions" (I'm loading requested  
> package:partitions)
> Errore: pacchetto 'partitions' non caricato (Error: package  
> "partitions" not loadad)
> Warning messages:
> 1: the package 'untb' has been created with R version 2.3.1
> 2: non c'è alcun pacchetto chiamato 'partitions' in: library(pkg,  
> character.only = TRUE, logical = TRUE, lib.loc = lib.loc) (there's  
> not a package named "partitions" in..)

This means that you need to also
install.packages('partitions')

You could have said
install.packages('untb', dependencies = TRUE),

and this would have done the job for you.
_
Professor Michael Kubovy
University of Virginia
Department of Psychology
USPS: P.O.Box 400400Charlottesville, VA 22904-4400
Parcels:Room 102Gilmer Hall
 McCormick RoadCharlottesville, VA 22903
Office:B011+1-434-982-4729
Lab:B019+1-434-982-4751
Fax:+1-434-982-4766
WWW:http://www.people.virginia.edu/~mk9y/

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

[R] Help with unlist

[Michael Kubovy]

Dear r-helpers,

I have a list whose elements are

 > str(durCut[[1]])
Ord.factor w/ 5 levels "vLow"<"low"<"med"<..: 3 2 5 2 2 2 4 4 3 5 ...

How do I unlist durCut into an ordered factor?

 > str(unlist(durCut))
int [1:3024] 3 2 5 2 2 2 4 4 3 5 ...

_____
Professor Michael Kubovy
University of Virginia
Department of Psychology
USPS: P.O.Box 400400Charlottesville, VA 22904-4400
Parcels:Room 102Gilmer Hall
 McCormick RoadCharlottesville, VA 22903
Office:B011+1-434-982-4729
Lab:B019+1-434-982-4751
Fax:+1-434-982-4766
WWW:http://www.people.virginia.edu/~mk9y/

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

[R] plot(all.effects)

[Michael Kubovy]

Dear r-helpers,

plot(all.effects()) uses trellis for its output, and often produces  
several effect plots in a layout. Can anyone tell me how to set the  
arguments to plot.effect which default to
row=1, col=1, nrow=1, ncol=1, more=FALSE
so that I can plot multiple effect plots on one layout?

I was hoping that a call like

last <- rep(c(TRUE, FALSE), c(3, 1))
lay <-matrix(c(1:4,), 2, 2, byrow = T)
for(i in 1:4) {
tmp.lm <- lm(y ~ a + b, data)
spl <- which(lay == i, arr.ind = TRUE)
rowN = spl[1]
colN = spl[2]
plot(all.effects(tmp.lm), row = rowN, col = colN, nrow = 2, ncol =  
2, more = last[i])
}

would work (my data are complicated, and have not had the time to set  
up a replicable example).

Here is what I ran, and what R told me:
for(i in 1:4) {
tmp.lm <- lm(duration ~ (pitch + vis + aud + display)^2 -  
aud:display - vis:aud -
pitch:vis, data = subset(txt, subj == unique(txt$subj)[i]))
spl <- which(lay == i, arr.ind = TRUE)
rowN = spl[1]
colN = spl[2]
plot(all.effects(tmp.lm), rug = F, ask = F,
main = paste('Subject = ',sep = '', i),
row = rown, col = coln, nrow = 2, ncol = 2, more = last[i])
  }
Error in plot.effect(x[[(i - 1) * cols + j]], row = i, col = j, nrow  
= rows,  :
formal argument "row" matched by multiple actual arguments

_____
Professor Michael Kubovy
University of Virginia
Department of Psychology
USPS: P.O.Box 400400Charlottesville, VA 22904-4400
Parcels:Room 102Gilmer Hall
 McCormick RoadCharlottesville, VA 22903
Office:B011+1-434-982-4729
Lab:B019+1-434-982-4751
Fax:+1-434-982-4766
WWW:http://www.people.virginia.edu/~mk9y/

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

[R] Constructing predictions from HPDinterval() after lmer()

[Michael Kubovy]

Dear r-helpers,

Following up on http://finzi.psych.upenn.edu/R/Rhelp02a/archive/ 
81159.html where Douglas Bates gives a helpful application of lmer()  
to data(sleepstudy, package = 'lme4'), I need a bit more help in  
order to plot the correct confidence intervals of a designed  
experiment such as:

 > data(ratdrink, package = 'faraway')

I follow the steps Douglas took in his example:

 > summary( rd.er <- lmer(wt ~ weeks*treat + (1 | subject), data =  
ratdrink) )

Linear mixed-effects model fit by REML
Formula: wt ~ weeks * treat + (1 | subject)
Data: ratdrink
AIC   BIC logLik MLdeviance REMLdeviance
962.4 982.8 -474.2  964.3948.4
Random effects:
Groups   NameVariance Std.Dev.
subject  (Intercept) 71.206   8.4384
Residual 51.222   7.1570
number of obs: 135, groups: subject, 27

Fixed effects:
   Estimate Std. Error t value
(Intercept)52.8800 3.1928   16.56
weeks  26.4800 0.7157   37.00
treatthiouracil 4.7800 4.51531.06
treatthyroxine -0.7943 4.9756   -0.16
weeks:treatthiouracil  -9.3700 1.0121   -9.26
weeks:treatthyroxine0.6629 1.11530.59

Correlation of Fixed Effects:
 (Intr) weeks  trtthr trtthy wks:trtthr
weeks   -0.448
treatthircl -0.707  0.317
treatthyrxn -0.642  0.288  0.454
wks:trtthrc  0.317 -0.707 -0.448 -0.203
wks:trtthyr  0.288 -0.642 -0.203 -0.448  0.454

 > rd.mc <- mcmcsamp(rd.er, 5)
 > library(coda)
 > HPDinterval(rd.mc)

lower  upper
(Intercept)46.420404  59.406398
weeks  25.070131  27.930363
treatthiouracil-4.420942  14.009291
treatthyroxine-10.758369   9.435761
weeks:treatthiouracil -11.404620  -7.337025
weeks:treatthyroxine   -1.603858   2.842704
log(sigma^2)3.683085   4.226153
log(sbjc.(In))  3.633853   4.955867
deviance  965.750351 980.825613
attr(,"Probability")
[1] 0.95


***To make sure my request is clear***
What I need is the analog of what is produced by all.effects() in the  
effects package:

 > summary(rd <- lm(wt ~ weeks*treat, data = ratdrink))

Call:
lm(formula = wt ~ weeks * treat, data = ratdrink)

Residuals:
 Min  1Q  Median  3Q Max
-23.514  -6.660   0.230   6.914  28.343

Coefficients:
   Estimate Std. Error t value Pr(>|t|)
(Intercept)52.8800 2.6547  19.919  < 2e-16 ***
weeks  26.4800 1.0838  24.433  < 2e-16 ***
treatthiouracil 4.7800 3.7544   1.2730.205
treatthyroxine -0.7943 4.1371  -0.1920.848
weeks:treatthiouracil  -9.3700 1.5327  -6.113 1.08e-08 ***
weeks:treatthyroxine0.6629 1.6890   0.3920.695
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Residual standard error: 10.84 on 129 degrees of freedom
Multiple R-Squared: 0.9121, Adjusted R-squared: 0.9087
F-statistic: 267.8 on 5 and 129 DF,  p-value: < 2.2e-16

 > rd.eff <- all.effects(rd)
 > rd.ci <- data.frame(y = rd.eff[[1]]$fit, Lower = rd.eff[[1]] 
$lower, Upper = rd.eff[[1]]$upper)

_
Professor Michael Kubovy
University of Virginia
Department of Psychology
USPS: P.O.Box 400400Charlottesville, VA 22904-4400
Parcels:Room 102Gilmer Hall
 McCormick RoadCharlottesville, VA 22903
Office:B011+1-434-982-4729
Lab:B019+1-434-982-4751
Fax:+1-434-982-4766
WWW:http://www.people.virginia.edu/~mk9y/

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

Re: [R] plotting 95% confidence bands on a simple linear regression model from lm()

[Michael Kubovy]

On Oct 20, 2006, at 2:46 PM, Jason Horn wrote:

> What's the best / simplest way to create 95% confidence bands for a
> model created with lm() that can be plotted around teh regression
> line?  I've looked everywhere for this - I guess I must be missing
> something.
>
> - Jason

library(effects)
?summary.effect
_________
Professor Michael Kubovy
University of Virginia
Department of Psychology
USPS: P.O.Box 400400Charlottesville, VA 22904-4400
Parcels:Room 102Gilmer Hall
 McCormick RoadCharlottesville, VA 22903
Office:B011+1-434-982-4729
Lab:B019+1-434-982-4751
Fax:+1-434-982-4766
WWW:http://www.people.virginia.edu/~mk9y/

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.