Re: [R] Lattice: regression lines within grouped xyplot panels

2008-08-09 Thread Rainer Hurling

On 09.08.2008 01:05 (UTC+1), Deepayan Sarkar wrote:

On Fri, Aug 8, 2008 at 2:38 PM, Rainer Hurling <[EMAIL PROTECTED]> wrote:

Dear community,

I am looking for a possibility to draw 'regression lines' instead of
'smooth' lines in grouped xyplots. The following code should give you a
small example of the data structure.


library(lattice)
data(Gcsemv, package = "mlmRev")

# Creates artificial grouping variable ...
Gcsemv$Groups <-
 ifelse(as.numeric(as.character(Gcsemv$school))>65000,
"Group1", "Group2")

xyplot(written ~ course | gender, data = Gcsemv,
  type = c("g", "p", "smooth"),
  groups = Groups,
  panel = function(x, y, ...) {
panel.xyplot(x, y, ...)
# Here I want to draw the regression lines
# panel.abline(x, y)
  },
  auto.key = list(space = 'right'))


Does this do what you want?:


Yes, exactly!


xyplot(written ~ course | gender, data = Gcsemv,
  type = c("g", "p", "r"),
  groups = Groups)

The problem with your approach is that the panel function you define
doesn't deal with groups. An easy workaround is to use
panel.superpose:


I knew that I had to look for a panel function dealing with groups, but 
I had no clue how to declare.


xyplot(written ~ course | gender, data = Gcsemv,
  type = c("g", "p"),
  groups = Groups,
  panel = panel.superpose,
  panel.groups = function(x, y, ...) {
panel.xyplot(x, y, ...)
panel.lmline(x, y, ...)
  },
  auto.key = list(space = 'right'))


Both, type("r") and panel.groups() are fine for my problem. And with 
panel.groups() I am able to write my own group functions, very nice. Now 
I can start analyzing ...  :-)



-Deepayan


Many thanks, also for your valuable book,
Rainer

__
R-help@r-project.org 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] import/export txt file

2008-08-09 Thread Dan Davison
On Fri, Aug 08, 2008 at 04:44:13PM -0700, Alessandro wrote:
> Hi All,
> 
>  
> 
> I have 2 questions:
> 
> 1.   Import: when I import my txt file (X,Y and Z) in R with "testground
> <- read.table(file="c:/work_LIDAR_USA/R_kriging/ground26841492694149.txt",
> header=T)", I lost the 4 number after the point (".").  does It possible add
> in the code the possibility to read the 4 numbers after the .
> 

I think the problem is simply that you have options()$digits set to 7
(the default). Read the 'digits' section in help(options) and try

options(digits=11)

>  
> 
> 2.   Does It possible to write a X, Y,  Z *txt file without the ID in R
> and sep"," for the rows?
> 

write.csv(your.data.frame, row.names=FALSE, quote=FALSE)

x <- read.table("path/to/your/file.txt", header=T)
x
#  X   Y   Z
# 1 26800.47 4149984 1543.39
# 2 26800.47 4149984 1543.39
options(digits=11)
x
#  X  Y   Z
# 1 26800.47 4149983.94 1543.39
# 2 26800.47 4149983.94 1543.39
> write.csv(x, row.names=FALSE, quote=FALSE)
# X,Y,Z
# 26800.47,4149983.94,1543.39
# 26800.47,4149983.94,1543.39

Dan

>  
> 
> Example:
> 
> Original data:
> 
> X Y Z
> 
> 26800.4700 4149983.9400 1543.3900
> 
> ... . ..
> 
>  
> 
> I wish to create a txt file (with "," sep):
> 
> X, Y, Z
> 
> 26800.4700, 4149983.9400, 1543.3900
> 
> ..., ., ..
> 
>  
> 
> Thanks (It's Friday night, sorry I am tired)
> 
>  
> 
> Ale
> 
> 
>   [[alternative HTML version deleted]]
> 
> __
> R-help@r-project.org 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@r-project.org 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] read.table question

2008-08-09 Thread Dan Davison
On Fri, Aug 08, 2008 at 07:27:13PM -0700, Alessandro wrote:
> Hi All.
> 
>  
> 
> I have a file txt with 3 columns (X, Y and Z).  every rows has 4 decimal
> place (i.e. x.). I use read.table to import the data in R, but with
> summary(), I don't see the decimal place after the dot. Is there any way for
> me to preserve the information?

I hope I've answered this in the first thread on the subject.

https://stat.ethz.ch/pipermail/r-help/2008-August/170422.html

Dan

p.s. People don't like it if you submit the same question twice.


> 
>  
> 
> testground <- read.table
> (file="c:/work_LIDAR_USA/R_kriging/ground26841492694149.txt", header=T)
> 
>  
> 
> thanks
> 
>  
> 
> Ale
> 
> 
>   [[alternative HTML version deleted]]
> 
> __
> R-help@r-project.org 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@r-project.org 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] levels values of cut()

2008-08-09 Thread baptiste auguie

Dear list,

 I have the following example, from which I am hoping to retrieve  
numeric values of the factor levels (that is, without the brackets):




x <- seq(1, 15, length=100)
y <- sin(x)

my.cuts <- cut(which(abs(y) < 1e-1), 3)
levels(my.cuts)


hist() does not suit me for this, as it does not necessarily respect  
the number of breaks.


getAnywhere hasn't got me very far: I cannot seem to find a readable  
code for the built-in cut function in the base library. I think  
getMethod should do it but I don't understand the arguments to pass.


Any pointers appreciated,

Thanks,

baptiste



_

Baptiste Auguié

School of Physics
University of Exeter
Stocker Road,
Exeter, Devon,
EX4 4QL, UK

Phone: +44 1392 264187

http://newton.ex.ac.uk/research/emag

__
R-help@r-project.org 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] levels values of cut()

2008-08-09 Thread Stephen Tucker
Not sure what you're looking for, but does this help?

Extending your code,

> library(gsubfn)
> t(strapply(levels(my.cuts),"([0-9.]+),([0-9.]+)",
+  function(...) as.numeric(c(...)),backref=-2,simplify=TRUE))
 [,1] [,2]
[1,] 15.9 38.3
[2,] 38.3 60.7
[3,] 60.7 83.1



- Original Message 
From: baptiste auguie <[EMAIL PROTECTED]>
To: r-help@r-project.org
Sent: Saturday, August 9, 2008 1:51:01 AM
Subject: [R] levels values of cut()

Dear list,

  I have the following example, from which I am hoping to retrieve  
numeric values of the factor levels (that is, without the brackets):

>
> x <- seq(1, 15, length=100)
> y <- sin(x)
>
> my.cuts <- cut(which(abs(y) < 1e-1), 3)
> levels(my.cuts)

hist() does not suit me for this, as it does not necessarily respect  
the number of breaks.

getAnywhere hasn't got me very far: I cannot seem to find a readable  
code for the built-in cut function in the base library. I think  
getMethod should do it but I don't understand the arguments to pass.

Any pointers appreciated,

Thanks,

baptiste



_

Baptiste Auguié

School of Physics
University of Exeter
Stocker Road,
Exeter, Devon,
EX4 4QL, UK

Phone: +44 1392 264187

http://newton.ex.ac.uk/research/emag

__
R-help@r-project.org 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@r-project.org 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] re cursive root finding

2008-08-09 Thread baptiste auguie

Hi all,

with a lot more fiddling around, I've come up with this solution:


x <- seq(1, 15, length=100)
y <- jitter(cos(x), a=0.2)

plot(x, y)

findZero <- function(x, y, guess = mean(range(x)), exclusion=  
5*diff(range(x))/length(x))

{
interval <- c(guess - exclusion, guess+exclusion)
print(interval)
   spl <- smooth.spline(x, y)
obj <- function(xx) {
z <- predict(spl, xx,  der=0)$y
sum((z)^2)
}
optimize(obj, interval)
}

# findZero(x, y) # test

findZeros <- function(x, y, ...){
spl <- smooth.spline(x, y)
y.spl <- predict(spl, x, der=0)$y
diff(diff(abs(y.spl), 1) > 0, 1)->test
guesses <- x[which(test==1)+1]
sapply(guesses, function(guess) findZero(x, y, guess, ...)$minimum)
}

zeros <- findZeros(x, y)
abline(v=zeros, col="red")


albeit not elegant or fast, nor very robust I suspect, it seems to  
work for me.


Best wishes

baptiste

On 8 Aug 2008, at 16:44, Hans W. Borchers wrote:



As your curve is defined by its points, I don't see any reason to
artificially apply functions such as 'uniroot' or 'optim' (being a
real overkill in this situation).

First smooth the curve with splines, Savitsky-Golay, or Whittacker
smoothing, etc., then loop through the sequence of points and compute
the gradient by hand, single-sided, two-sided, or both.

At the same time, mark those indices where the gradient is zero or
changes its sign; these will be the solutions you looked for.

With your example, I immediate got as maxima or minima:

   x1 = 1.626126
   x2 = 4.743243
   x3 = 7.851351

//  Hans Werner Borchers


Any comments? Maybe the problem was not clear or looked too specific.
I'll add a more "bare bones" example, if only to simulate discussion:


x <- seq(1,10,length=1000)
set.seed(11) # so that you get the same randomness
y <- jitter(sin(x),a = 0.2)
values <- data.frame(x= x,  y = y)

findZero <- function (guess, neighbors, values)
{

   smooth <- smooth.spline(values$x, values$y)

   obj <- function(x) {
   (predict(smooth, x)$y) ^2
   }
   minimum <- which.min(abs(values$x - guess))
   optimize(obj, interval = c(values$x[minimum - neighbors],
   values$x[minimum + 
neighbors]))  # uniroot could be used
instead i suppose

}

test <- findZero(guess = 6 ,  neigh = 50, values = values)
plot(x,y)
abline(h=0)
abline(v=test$minimum, col="red")


Now, I would like to find all (N=)3 roots, without a priori knowledge
of their location in the interval. I considered several approaches:

1) find all the numerical values of the initial data that are close to
zero with a given threshold. Group these values in N sets using cut()
and hist() maybe? I could never get this to work, the factors given by
cut confuse me (i couldn't extract their value). Then, apply the
function given above with the guess given by the center of mass of the
different groups of zeros.

2) apply findZero once, store the result, then add something big
(1e10) to the neighboring points and look for a zero again and repeat
the procedure until N are found. This did not work, I assume because
it does not perturb the minimization problem in the way I want.

3) once a zero is found, look for zeros on both sides, etc... this
quickly makes a complicated decision tree when the number of zeros
grows and I could not find a clean way to implement it.

Any thoughts welcome! I feel like I've overlooked an obvious trick.

Many thanks,

baptiste


On 7 Aug 2008, at 11:49, baptiste auguie wrote:


Dear list,


I've had this problem for a while and I'm looking for a more general
and robust technique than I've been able to imagine myself. I need
to find N (typically N= 3 to 5) zeros in a function that is not a
polynomial in a specified interval.

The code below illustrates this, by creating a noisy curve with
three peaks of different position, magnitude, width and asymmetry:


x <- seq(1, 10, length=500)
exampleFunction <- function(x){ # create some data with peaks of
different scaling and widths + noise
fano <- function (p, x)
{
y0 <- p[1]
A1 <- abs(p[2])
w1 <- abs(p[3])
x01 <- p[4]
q <- p[5]
B1 <- (2 * A1/pi) * ((q * w1 + x - x01)^2/(4 * (x - x01)^2 +
w1^2))
y0 + B1
}
p1 <- c(0.1, 1, 1, 5, 1)
p2 <- c(0.5, 0.7, 0.2, 4, 1)
p3 <- c(0, 0.5, 3, 1.2, 1)
y <- fano(p1, x) + fano(p2, x) + fano(p3, x)
jitter(y, a=0.005*max(y))
}

y <- exampleFunction(x)

sample.df <- data.frame(x = x, y = y)

with(sample.df, plot(x, y, t="l")) # there are three peaks, located
around x=2, 4 ,5
y.spl <- smooth.spline(x, y) # smooth the noise
lines(y.spl, col="red")



I wish to obtain the zeros of the first and second derivatives of
the smoothed spline y.spl. I can use uniroot() or optim() to find
one

Re: [R] RPro

2008-08-09 Thread Ajay ohri
I think a bit of standardization and commercial exposure would help R, given
the fact it is much much better than any other expensive stats package yet
has no gotten the market share it deserves because everyone is either
building it or using it, hence no one has actually taken the effort to
promote it.

Also R package developers are overwhelmingly for  open source so they would
be able to resist any "devious" commercial angle.

Ajay

www.decisionstats.com



On Fri, Aug 8, 2008 at 11:55 PM, Kenn Konstabel <[EMAIL PROTECTED]> wrote:

> There's more to this trend: SPSS and Statistica now advertise "R language
> support" :
>
> http://www.statsoft.com/industries/Rlanguage.htm
> http://www.spss.com/spssdirections/na/sessions.cfm?sessionType=2
>
>
>
> Kenn Konstabel
>
> On Fri, Aug 8, 2008 at 5:30 PM, Marc Schwartz <[EMAIL PROTECTED]
> >wrote:
>
> > on 08/08/2008 09:13 AM [EMAIL PROTECTED] wrote:
> >
> >> I recently came across a flyer from REvolution Computing, and I wanted
> to
> >> ask if this is R going private?
> >> Tony.
> >>
> >
> > It is one of at least two commercial offerings of R that are being
> > developed/released. The other, that I know of, is RStat:
> >
> >  http://www.random-technologies-llc.com/
> >
> > As was recently discussed in a lengthy thread, this is entirely permitted
> > under the terms of the GPL.
> >
> > HTH,
> >
> > Marc Schwartz
> >
> >
> > __
> > R-help@r-project.org 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.
> >
>
> [[alternative HTML version deleted]]
>
> __
> R-help@r-project.org 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.
>

[[alternative HTML version deleted]]

__
R-help@r-project.org 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] levels values of cut()

2008-08-09 Thread Prof Brian Ripley

On Sat, 9 Aug 2008, baptiste auguie wrote:


Dear list,

I have the following example, from which I am hoping to retrieve numeric 
values of the factor levels (that is, without the brackets):




x <- seq(1, 15, length=100)
y <- sin(x)

my.cuts <- cut(which(abs(y) < 1e-1), 3)
levels(my.cuts)


hist() does not suit me for this, as it does not necessarily respect the 
number of breaks.


getAnywhere hasn't got me very far: I cannot seem to find a readable code for 
the built-in cut function in the base library. I think getMethod should do it 
but I don't understand the arguments to pass.


Not getMethod (that's for S4 methods).  Just type cut.default at the R 
prompt.


However, try

example(cut)
foo <- levels(cut(aaa, 3))
lims <- matrix(nrow=length(foo), ncol=2)
lims[,1] <- as.numeric( sub("\\((.+),.*", "\\1", foo) )
lims[,2] <- as.numeric( sub("[^,]*,([^]]*)\\]", "\\1", foo) )

--
Brian D. Ripley,  [EMAIL PROTECTED]
Professor of Applied Statistics,  http://www.stats.ox.ac.uk/~ripley/
University of Oxford, Tel:  +44 1865 272861 (self)
1 South Parks Road, +44 1865 272866 (PA)
Oxford OX1 3TG, UKFax:  +44 1865 272595

__
R-help@r-project.org 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] effective matrix subset

2008-08-09 Thread jgarcia
Hi;
If we have a matrix A, and a vector X, where length(X)=nrow(A), and X
contains a wanted column for each row in A, in row ascending order. How
would be the most effective way to extract the desired vector V (with
length(V)=nrow(A))?

Wishes,
Javier

__
R-help@r-project.org 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] effective matrix subset

2008-08-09 Thread Marc Schwartz

on 08/09/2008 06:01 AM [EMAIL PROTECTED] wrote:

Hi;
If we have a matrix A, and a vector X, where length(X)=nrow(A), and X
contains a wanted column for each row in A, in row ascending order. How
would be the most effective way to extract the desired vector V (with
length(V)=nrow(A))?



A <- matrix(1:20, 4, 5)

> A
 [,1] [,2] [,3] [,4] [,5]
[1,]159   13   17
[2,]26   10   14   18
[3,]37   11   15   19
[4,]48   12   16   20


# Create an arbitrary set of indices, one for each row in A
X <- c(2, 5, 1, 4)

> X
[1] 2 5 1 4


Presumably you want:

V <- c(A[1, 2], A[2, 5], A[3, 1], A[4, 4])

> V
[1]  5 18  3 16


If so, then:

> sapply(seq(nrow(A)), function(i) A[i, X[i]])
[1]  5 18  3 16


Is that what you were looking for?


BTW, see ?diag for a special case:

> diag(A)
[1]  1  6 11 16


HTH,

Marc Schwartz

__
R-help@r-project.org 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] effective matrix subset

2008-08-09 Thread Dan Davison
On Sat, Aug 09, 2008 at 06:29:59AM -0500, Marc Schwartz wrote:
> on 08/09/2008 06:01 AM [EMAIL PROTECTED] wrote:
>> Hi;
>> If we have a matrix A, and a vector X, where length(X)=nrow(A), and X
>> contains a wanted column for each row in A, in row ascending order. How
>> would be the most effective way to extract the desired vector V (with
>> length(V)=nrow(A))?
>
>
> A <- matrix(1:20, 4, 5)
>
> > A
>  [,1] [,2] [,3] [,4] [,5]
> [1,]159   13   17
> [2,]26   10   14   18
> [3,]37   11   15   19
> [4,]48   12   16   20
>
>
> # Create an arbitrary set of indices, one for each row in A
> X <- c(2, 5, 1, 4)
>
> > X
> [1] 2 5 1 4
>
>
> Presumably you want:
>
> V <- c(A[1, 2], A[2, 5], A[3, 1], A[4, 4])
>
> > V
> [1]  5 18  3 16
>
>
> If so, then:
>
> > sapply(seq(nrow(A)), function(i) A[i, X[i]])
> [1]  5 18  3 16

Or

> A[cbind(seq(nrow(A)), X)]
[1]  5 18  3 16

Dan

>
>
> Is that what you were looking for?
>
>
> BTW, see ?diag for a special case:
>
> > diag(A)
> [1]  1  6 11 16
>
>
> HTH,
>
> Marc Schwartz
>
> __
> R-help@r-project.org 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@r-project.org 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] [lme4]Coef output with binomial lmer

2008-08-09 Thread Mark Difford

Hi Tom,

>> 1|ass%in%pop%in%fam

This is "non-standard," but as you have found, it works. The correct
translation is in fact

1|fam/pop/ass

and not 1|ass/pop/fam as suggested by Harold Doran. Dropping %,
ass%in%pop%in%fam reads [means] as: nest ass in pop [= pop/ass], and then
nest this in fam == fam/pop/ass

HTH, Mark.


T.C. Cameron wrote:
> 
> Dear R users  
>  
> I have built the following model
>  
> m1<-lmer(y~harn+foodn+(1|ass%in%pop%in%fam),family = "quasibinomial")
>  
> where y<-cbind(alive,dead)
>  
> where harn and foodn are categorical factors and the random effect is a
> nested term to represent experimental structure
> e.g.  Day/Block/Replicate
> ass= 5 level factor, pop= 2 populations per treatment factor in each
> assay, 7 reps per population
>  
> The model can be family = quasibinomial or binomial
>  
> My complete lack of understanding is in retrieving the coefficients for
> the fixed effects to back-transform the effects of my factors on
> proportional survival
>  
> I get the following output:
>> coef(m1)
> $`ass %in% pop %in% fam`
>   (Intercept)  harn1 harn2   foodn2
> FALSE   1.0322375 -0.1939521 0.0310434 0.810084
> TRUE0.5997679 -0.1939521 0.0310434 0.810084
>  
> Where FALSE and TRUE refer to some attribute of the random effect 
>  
> My hunch is that it refers to the Coefficients with (=TRUE) and without
> (=FALSE) the random effects?
>  
> Any help appreciated
>  
> 
> 
> 
> Dr Tom C Cameron
> Genetics, Ecology and Evolution
> IICB, University of Leeds
> Leeds, UK
> Office: +44 (0)113 343 2837
> Lab:+44 (0)113 343 2854
> Fax:+44 (0)113 343 2835
> 
> 
> Email: [EMAIL PROTECTED]
> Webpage: click here
>  
> 
>  
> 
>   [[alternative HTML version deleted]]
> 
> __
> R-help@r-project.org 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.
> 
> 

-- 
View this message in context: 
http://www.nabble.com/-lme4-Coef-output-with-binomial-lmer-tp18894407p18904468.html
Sent from the R help mailing list archive at Nabble.com.

__
R-help@r-project.org 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] Tick marks that correspond with bars on barplot

2008-08-09 Thread Jim Lemon
On Fri, 2008-08-08 at 16:04 -0300, Megan J Bellamy wrote:
> Hello all,
> 
> I have created a barplot that shows change in hardwood/softwood density from 
> 1965 to 2005 in 5 year periods (1965,1970, etc). I would like to have an 
> X-axis where the labels for the years line up after every two bars in the 
> plot (there is one bar for hardwood, and another for softwood). Below is my 
> script:
> 
> density<-read.table("F:\\Megan\\Vtest.csv", header=TRUE, sep=",")
> attach (density)
> barplot(DENSITY,YEAR, col=c("blue", "green", "green", "blue", "blue", 
> "green", "blue", "green", "green", "blue", "green", "blue", "blue", "green" , 
> "green", "blue", "blue", "green"))
> legend(1,85,c("Softwood", "Hardwood"), fill=c("blue", "green"), bty="n")
> title(ylab="Density of Softwood and Hardwood by percent (%)")
> title(xlab="Year of measurement")
> title(main="Change in Softwood and Hardwood Density between 1965-2005 for PSP 
> #80")
> axis(2, at=NULL)
> 
> I have tried using the tck, axis.lty functions with no luck and have also 
> tried
> axis(1, at=YEAR) # but only the first year (1965) comes up. 
> 
Try this:

xpos<-barplot(...)
...
axis(1,at=xpos,labels=seq(1965,2005,by=5))

> Attached is the csv file. Any help with this would be greatly appreciated. 
> Many thanks in advance,
> 
Unfortunately, the CSV file didn't make it, but the above works with
simulated data.

Jim

__
R-help@r-project.org 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] levels values of cut()

2008-08-09 Thread baptiste auguie
Thank you all for the precious tips. For memory I've made the  
following wrapper function for this. I wonder whether a short note on  
these regular expressions could be useful on the help page of cut().




cutIntervals <- function(x, ...){
dotArgs <- unlist(c(...))
	if( any(names(dotArgs) == "labels")) stop("labels cannot be  
specified,  use cut instead")


cut.fact <- levels(cut(x,labels=NULL, ...))
# tip from Brian Ripley
lims <- matrix(nrow=length(cut.fact), ncol=2)
lims[,1] <- as.numeric( sub("\\((.+),.*", "\\1", cut.fact) )
lims[,2] <- as.numeric( sub("[^,]*,([^]]*)\\]", "\\1", cut.fact) )
# alternatively (Stephen Tucker)
 # library(gsubfn)
 # lims <- t(strapply(cut.fact,"([0-9.]+),([0-9.]+)",
 #  function(...) 
as.numeric(c(...)),backref=-2,simplify=TRUE))
lims
}

cutIntervals(1:5, 3)



Many thanks,

baptiste

On 9 Aug 2008, at 11:12, Prof Brian Ripley wrote:


On Sat, 9 Aug 2008, baptiste auguie wrote:


Dear list,

I have the following example, from which I am hoping to retrieve  
numeric values of the factor levels (that is, without the brackets):



x <- seq(1, 15, length=100)
y <- sin(x)
my.cuts <- cut(which(abs(y) < 1e-1), 3)
levels(my.cuts)


hist() does not suit me for this, as it does not necessarily  
respect the number of breaks.


getAnywhere hasn't got me very far: I cannot seem to find a  
readable code for the built-in cut function in the base library. I  
think getMethod should do it but I don't understand the arguments  
to pass.


Not getMethod (that's for S4 methods).  Just type cut.default at the  
R prompt.


However, try

example(cut)
foo <- levels(cut(aaa, 3))
lims <- matrix(nrow=length(foo), ncol=2)
lims[,1] <- as.numeric( sub("\\((.+),.*", "\\1", foo) )
lims[,2] <- as.numeric( sub("[^,]*,([^]]*)\\]", "\\1", foo) )

--
Brian D. Ripley,  [EMAIL PROTECTED]
Professor of Applied Statistics,  http://www.stats.ox.ac.uk/~ripley/
University of Oxford, Tel:  +44 1865 272861 (self)
1 South Parks Road, +44 1865 272866 (PA)
Oxford OX1 3TG, UKFax:  +44 1865 272595


_

Baptiste Auguié

School of Physics
University of Exeter
Stocker Road,
Exeter, Devon,
EX4 4QL, UK

Phone: +44 1392 264187

http://newton.ex.ac.uk/research/emag

__
R-help@r-project.org 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] effective matrix subset

2008-08-09 Thread Marc Schwartz

on 08/09/2008 06:52 AM Dan Davison wrote:

On Sat, Aug 09, 2008 at 06:29:59AM -0500, Marc Schwartz wrote:

on 08/09/2008 06:01 AM [EMAIL PROTECTED] wrote:

Hi;
If we have a matrix A, and a vector X, where length(X)=nrow(A), and X
contains a wanted column for each row in A, in row ascending order. How
would be the most effective way to extract the desired vector V (with
length(V)=nrow(A))?


A <- matrix(1:20, 4, 5)


A

 [,1] [,2] [,3] [,4] [,5]
[1,]159   13   17
[2,]26   10   14   18
[3,]37   11   15   19
[4,]48   12   16   20


# Create an arbitrary set of indices, one for each row in A
X <- c(2, 5, 1, 4)


X

[1] 2 5 1 4


Presumably you want:

V <- c(A[1, 2], A[2, 5], A[3, 1], A[4, 4])


V

[1]  5 18  3 16


If so, then:


sapply(seq(nrow(A)), function(i) A[i, X[i]])

[1]  5 18  3 16


Or


A[cbind(seq(nrow(A)), X)]

[1]  5 18  3 16

Dan


Better (and faster) solution Dan.

I can't blame the lack of coffee on missing that one this morning. I 
have had a full pot already over the past 6 hours, working on shifting 
my internal clock and getting ready to begin my journey to Dortmund 
later tonight...


Safe travels to all who are going.

Marc

__
R-help@r-project.org 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] Color of box frame in Legend (Was: Matrix barplot)

2008-08-09 Thread Andreas Tille

On Mon, 28 Jul 2008, Nutter, Benjamin wrote:


Try sourcing in the 'new.legend' function below. It's the legend
function with a new argument called 'box.col'.  The argument will change
the color of the box surrounding the legend.  If I understand what it is
you are looking for, this should work.  Also, I didn't see a way to
change the axis bar in your code, so I suppressed the axis in the call
to barplot, and manually replaced the axis using the axis function.  I
hope this works for you.


Sorry for replying late - I was busy with traveling to DebConf in Argentina
and doing important stuff here.  No, your sulution does not make the situation
better - it actually makes it worse because also the lines around the graphs
become black instead of giving me some influence on the lines around the
boxes in the legend.

Kind regards

   Andreas.

--
http://fam-tille.de

__
R-help@r-project.org 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] levels values of cut()

2008-08-09 Thread Prof Brian Ripley

On Sat, 9 Aug 2008, baptiste auguie wrote:

Thank you all for the precious tips. For memory I've made the following 
wrapper function for this. I wonder whether a short note on these regular 
expressions could be useful on the help page of cut().


Already there in R-devel 





cutIntervals <- function(x, ...){
dotArgs <- unlist(c(...))
	if( any(names(dotArgs) == "labels")) stop("labels cannot be 
specified,  use cut instead")


cut.fact <- levels(cut(x,labels=NULL, ...))
# tip from Brian Ripley
lims <- matrix(nrow=length(cut.fact), ncol=2)
lims[,1] <- as.numeric( sub("\\((.+),.*", "\\1", cut.fact) )
lims[,2] <- as.numeric( sub("[^,]*,([^]]*)\\]", "\\1", cut.fact) )
# alternatively (Stephen Tucker)
 # library(gsubfn)
 # lims <- t(strapply(cut.fact,"([0-9.]+),([0-9.]+)",
	 # 		function(...) 
as.numeric(c(...)),backref=-2,simplify=TRUE))

lims
}

cutIntervals(1:5, 3)



Many thanks,

baptiste

On 9 Aug 2008, at 11:12, Prof Brian Ripley wrote:


On Sat, 9 Aug 2008, baptiste auguie wrote:


Dear list,

I have the following example, from which I am hoping to retrieve numeric 
values of the factor levels (that is, without the brackets):



x <- seq(1, 15, length=100)
y <- sin(x)
my.cuts <- cut(which(abs(y) < 1e-1), 3)
levels(my.cuts)


hist() does not suit me for this, as it does not necessarily respect the 
number of breaks.


getAnywhere hasn't got me very far: I cannot seem to find a readable code 
for the built-in cut function in the base library. I think getMethod 
should do it but I don't understand the arguments to pass.


Not getMethod (that's for S4 methods).  Just type cut.default at the R 
prompt.


However, try

example(cut)
foo <- levels(cut(aaa, 3))
lims <- matrix(nrow=length(foo), ncol=2)
lims[,1] <- as.numeric( sub("\\((.+),.*", "\\1", foo) )
lims[,2] <- as.numeric( sub("[^,]*,([^]]*)\\]", "\\1", foo) )

--
Brian D. Ripley,  [EMAIL PROTECTED]
Professor of Applied Statistics,  http://www.stats.ox.ac.uk/~ripley/
University of Oxford, Tel:  +44 1865 272861 (self)
1 South Parks Road, +44 1865 272866 (PA)
Oxford OX1 3TG, UKFax:  +44 1865 272595


_

Baptiste Auguié

School of Physics
University of Exeter
Stocker Road,
Exeter, Devon,
EX4 4QL, UK

Phone: +44 1392 264187

http://newton.ex.ac.uk/research/emag

__
R-help@r-project.org 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.


--
Brian D. Ripley,  [EMAIL PROTECTED]
Professor of Applied Statistics,  http://www.stats.ox.ac.uk/~ripley/
University of Oxford, Tel:  +44 1865 272861 (self)
1 South Parks Road, +44 1865 272866 (PA)
Oxford OX1 3TG, UKFax:  +44 1865 272595__
R-help@r-project.org 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] Reshape set operations?

2008-08-09 Thread rkevinburton
I have mange to use the library reshape to give me data structures that I want. 
Specifically:

m2008 <- melt(t2008, id.var=c("DayOfYear","Category","SubCategory","Sku"), 
measure.var=c("Quantity"))
 m2007 <- melt(t2007, id.var=c("DayOfYear","Category","SubCategory","Sku"), 
measure.var=c("Quantity"))

 r2008 <- cast(m2008, DayOfYear ~ variable | Sku, sum)
 r2007 <- cast(m2007, DayOfYear ~ variable | Sku, sum)

Now I would like to union the two lists. So I start out with an empty master 
list that will contain (when I am done) the merge (union) of r2008 and r2007. 
By "union" I mean that if the Sku exists in r2007 and r2008 I would like to 
create a new data frame that has the lists for DayOfYear and Quantity "merged" 
and append it to 'master'. If the Sku is not common to both objects then just 
copy or append to the 'master'.

Is it possible to come up with an expression and aggregate function that will 
do this? Is this better handle with aggregate functions or lapply? Being new to 
'R' it is hard for me to tell.

Thank you.

Kevin

__
R-help@r-project.org 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] RPro

2008-08-09 Thread Charles C. Berry

On Fri, 8 Aug 2008, Kenn Konstabel wrote:


There's more to this trend: SPSS and Statistica now advertise "R language
support" :

http://www.statsoft.com/industries/Rlanguage.htm
http://www.spss.com/spssdirections/na/sessions.cfm?sessionType=2



If you can't beat R, join R.

---

Should someone start a pool on when SAS will offer 'R language support' ??

;-)

Chuck






Kenn Konstabel

On Fri, Aug 8, 2008 at 5:30 PM, Marc Schwartz <[EMAIL PROTECTED]>wrote:


on 08/08/2008 09:13 AM [EMAIL PROTECTED] wrote:


I recently came across a flyer from REvolution Computing, and I wanted to
ask if this is R going private?
Tony.



It is one of at least two commercial offerings of R that are being
developed/released. The other, that I know of, is RStat:

 http://www.random-technologies-llc.com/

As was recently discussed in a lengthy thread, this is entirely permitted
under the terms of the GPL.

HTH,

Marc Schwartz


__
R-help@r-project.org 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.



[[alternative HTML version deleted]]

__
R-help@r-project.org 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.



Charles C. Berry(858) 534-2098
Dept of Family/Preventive Medicine
E mailto:[EMAIL PROTECTED]  UC San Diego
http://famprevmed.ucsd.edu/faculty/cberry/  La Jolla, San Diego 92093-0901

__
R-help@r-project.org 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] RPro

2008-08-09 Thread Marc Schwartz

on 08/09/2008 12:13 PM Charles C. Berry wrote:

On Fri, 8 Aug 2008, Kenn Konstabel wrote:


There's more to this trend: SPSS and Statistica now advertise "R language
support" :

http://www.statsoft.com/industries/Rlanguage.htm
http://www.spss.com/spssdirections/na/sessions.cfm?sessionType=2



If you can't beat R, join R.


"Resistance is futile. You will be assimilated."

;-)



Should someone start a pool on when SAS will offer 'R language support' ??


It is already in place, albeit not from Cary:

  http://www.minequest.com/

See the link for "A Bridge to R for SAS Users".

That product might put an interesting spin on:

> fortune(28)

Jim Gustafsson: I would like to put my SAS-code into R. Could I do
that, if yes, how?
Frank Harrell: Just reverse the procedure you use when you put R code
into SAS. ;)
   -- Jim Gustafsson and Frank Harrell
  R-help (February 2004)


Their WPS product has been getting a fair amount of attention, at least 
for those needing mostly "Base" and not "Stat" functionality at 
substantially less cost than the SAS commercial license pricing.


Regards,

Marc

__
R-help@r-project.org 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] RPro

2008-08-09 Thread Marc Schwartz

on 08/09/2008 12:45 PM Marc Schwartz wrote:

on 08/09/2008 12:13 PM Charles C. Berry wrote:

On Fri, 8 Aug 2008, Kenn Konstabel wrote:

There's more to this trend: SPSS and Statistica now advertise "R 
language

support" :

http://www.statsoft.com/industries/Rlanguage.htm
http://www.spss.com/spssdirections/na/sessions.cfm?sessionType=2



If you can't beat R, join R.


"Resistance is futile. You will be assimilated."

;-)


Should someone start a pool on when SAS will offer 'R language 
support' ??


It is already in place, albeit not from Cary:

  http://www.minequest.com/

See the link for "A Bridge to R for SAS Users".

That product might put an interesting spin on:

 > fortune(28)

Jim Gustafsson: I would like to put my SAS-code into R. Could I do
that, if yes, how?
Frank Harrell: Just reverse the procedure you use when you put R code
into SAS. ;)
   -- Jim Gustafsson and Frank Harrell
  R-help (February 2004)


Their WPS product has been getting a fair amount of attention, at least 
for those needing mostly "Base" and not "Stat" functionality at 
substantially less cost than the SAS commercial license pricing.


Correction. WPS is actually a product of "Word Programming", which is at:

  http://www.teamwpc.co.uk/products/wps

though at the moment, I am having trouble getting to their web page. 
There is a Wikipedia page here:


  http://en.wikipedia.org/wiki/World_Programming_System

Minequest is providing SAS to WPS related consulting services.

Marc

__
R-help@r-project.org 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] Box.test degrees of freedom

2008-08-09 Thread David Stoffer

I believe tsdiag() uses the correct degrees of freedom in applying Box.test,
but the graphic shows "lag" on the horizontal axis when it should display
"degrees of freedom".   



raf.rossignol wrote:
> 
> Hello,
> 
> Prof Brian Ripley wrote:
>> 
>> I think you are referring to its application to the residuals of an 
>> ARMA(p, q) fit, and that is not what Box.test says it does.
>> 
>> It is very easy to edit the code if you want to use a different degrees
>> of
>> freedom.
>> 
> I am also new to R, but it seems to me that there is still something
> confusing, not in Box.test but in tsdiag.Arima
> Indeed, the help of tsdiag says "The methods for 'arima' and 'StructTS'
> [...] use the Ljung-Box version of the portmanteau test." 
> So we could expect the degrees of freedom 'h-p-q' to be used, but a look
> at tsdiag.Arima shows it uses Box.test at lags h=1:gof.lag, with degrees
> of freedom equal to h, and  not h-p-q. Do you think this is a mistake in
> tsdiag.Arima or is there some experimental (or theoretical) reason
> supporting this choice ?
> Best regards
> 
> 
> 


-
The power of accurate observation is commonly called cynicism 
by those who have not got it.  George Bernard Shaw
-- 
View this message in context: 
http://www.nabble.com/Box.test-degrees-of-freedom-tp17277964p18907921.html
Sent from the R help mailing list archive at Nabble.com.

__
R-help@r-project.org 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 interpret t.test output

2008-08-09 Thread Felipe Carrillo
# Hi all:
  #I got a vector with fish lengths(mm)
 # Can someone help me interpret the output of
 # a t.test in plain english?
 #  Based on the t.test below I can say that
 # I reject the null hypothesis because
 # the p-value is smaller than the the significance
 # level(alpha=0.05). What else can I conclude here?
Ho = 36 mm
Ha <> 36 mm
fishlength <- c(35,32,37,39,42,45,37,36,35,34,40,42,41,50)
t.test(fishlength,mu=36)
   One Sample t-test

data:  fishlength
t = 2.27, df = 13, p-value = 0.04087
alternative hypothesis: true mean is not equal to 36
95 percent confidence interval:
 36.14141 41.71573
sample estimates:
mean of x
 38.92857

# I also would like to know how to interpret the output
# when mu=0? Notice that the p-value from the t.test above
#is different from this t.test. Are we trying to reject the null 
# hypothesis here too?

t.test(fishlength)
 One Sample t-test

data:  fishlength
t = 30.1741, df = 13, p-value = 2.017e-13
alternative hypothesis: true mean is not equal to 0
95 percent confidence interval:
 36.14141 41.71573
sample estimates:
mean of x
 38.92857

Thanks in advance for your help.


Felipe D. Carrillo
Supervisory Fishery Biologist  
Department of the Interior  
US Fish & Wildlife Service  
California, USA

__
R-help@r-project.org 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] Multivariate Regression with Weights

2008-08-09 Thread Arne Henningsen
> -Original Message-
> From: [EMAIL PROTECTED] [mailto:[EMAIL PROTECTED]
> Sent: Monday, August 04, 2008 5:15 PM
> To: Zhang Yanwei - Princeton-MRAm
> Subject: RE: [R] Multivariate Regression with Weights
>
>   the systemfit package can do that and the documentation is quite nice
> also.

The systemfit() function of the systemfit package can use different weights 
for the different *equations* but NOT different weights for different 
*observations* (if argument "method" is "SUR", "WLS", "3SLS" or "W2SLS").
I guess that you can divide the first equation by x_1 and the second equation 
by x_2 to obtain
   y_1^* = 1 + x_2^* + e_1^*
   y_2^* = 1 + x_1^* + e_2^*
with 
   y_1^* = y_1 / x_1
   x_2^* = x_2 / x_1
   e_1^* = e_1 / x_1
   y_2^* = y_2 / x_2
   x_1^* = x_1 / x_2
   e_2^* = e_2 / x_2

In this case, you should have homoscedastic error terms. Hence, you can use 
systemfit to estimate the system of the modified equations.

> On Mon, Aug 4, 2008 at  4:39 PM, Zhang Yanwei - Princeton-MRAm wrote:
> > Hi all,
> >   I'd like to fit a multivariate regression with the variance of the
> > error term porportional to the predictors, like the WLS in the
> > univariate case.
> >   y_1~x_1+x_2
> >   y_2~x_1+x_2
> >   var(y_1)=x_1*sigma_1^2

I guess that you mean var(e_1). Right?

> >   var(y_2)=x_2*sigma_2^2
> >   cov(y_1,y_2)=sqrt(x_1*x_2)*sigma_12^2
> >
> >  How can I specify this in R? Is there a corresponding function to the
> > univariate specification lm(y~x,weights=x)?? Thanks.

Unfortunately, systemfit has not (yet) a "weights" argument. (If you add this 
feature, I would be happy if you send me the code that I can include this 
feature in the official version.)

Best wishes,
Arne

-- 
Arne Henningsen
[EMAIL PROTECTED]
http://www.arne-henningsen.name


pgpDWkApw1m0r.pgp
Description: PGP signature
__
R-help@r-project.org 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 interpret t.test output

2008-08-09 Thread Ted Harding
On 09-Aug-08 20:31:33, Felipe Carrillo wrote:
># Hi all:
>  #I got a vector with fish lengths(mm)
>  # Can someone help me interpret the output of
>  # a t.test in plain english?
>  #  Based on the t.test below I can say that
>  # I reject the null hypothesis because
>  # the p-value is smaller than the the significance
>  # level(alpha=0.05). What else can I conclude here?
> Ho = 36 mm
> Ha <> 36 mm
> fishlength <- c(35,32,37,39,42,45,37,36,35,34,40,42,41,50)
> t.test(fishlength,mu=36)
>One Sample t-test
> 
> data:  fishlength
> t = 2.27, df = 13, p-value = 0.04087
> alternative hypothesis: true mean is not equal to 36
> 95 percent confidence interval:
>  36.14141 41.71573
> sample estimates:
> mean of x
>  38.92857

The standard interpretation of the above is that, on the hypothesis
that the fish were sampled from a population in which the mean length
was 36 (and the distribution of length is sufficiently close to
Normal for the t-test to be adequately applicable), then the value
of t differs from 0 by an amount which has probability less than
0.05 of being attained if that hypothesis were true. Therefore the
result is evidence of some strength that the hypotesis is not true,
and the mean length of fish in the population is different from 36.

However, it is on the borderline of the mildest criterion of
"significant difference" usually adopted, namely 0.05. A more
stringent criterion could be, for instance, 0.01; and then you
would not be "rejecting this null hypothesis".

As for the following:

> # I also would like to know how to interpret the output
> # when mu=0? Notice that the p-value from the t.test above
> #is different from this t.test. Are we trying to reject the null 
> # hypothesis here too?
> 
> t.test(fishlength)
>  One Sample t-test
> 
> data:  fishlength
> t = 30.1741, df = 13, p-value = 2.017e-13
> alternative hypothesis: true mean is not equal to 0
> 95 percent confidence interval:
>  36.14141 41.71573
> sample estimates:
> mean of x
>  38.92857
> 
> Thanks in advance for your help.

In terms of interpreting a statistical test, using your data,
of the hypothesis that the mean length in the population is 0,
the P-value of 0.2017 is very strong evidence indeed
that the mean is not 0.

However, I do not know why you are asking the question. No test
is needed. The length of any living fish, even while it is still
in the egg, is greater than 0; and whatever population you have
taken your sample from will have a mean length which is greater
than 0.

That is not to say that the result of a t-test on any sample
will necessarily give a significant result. You could have a
small catch with lengths, say,

   fishlengths <- c(2,4,9,20,50)
   t.test(fishlengths,mu=0)

# One Sample t-test
# data:  fishlengths 
# t = 1.9273, df = 4, p-value = 0.1262
# alternative hypothesis: true mean is not equal to 0 
# 95 percent confidence interval:
#  -7.489442 41.489442 
# sample estimates:
# mean of x 
#17 

And all you can conlude from that is that the sample, *in itself*,
does not carry sufficient information to confirm what you know
is true (i.e. mu > 0). Even the one-sided test of mu=0 with alternative
alt="greater" does not give a result significant at 5%:

 t.test(fishlengths,mu=0,alt="greater")

# One Sample t-test
# data:  fishlengths 
# t = 1.9273, df = 4, p-value = 0.0631
# alternative hypothesis: true mean is greater than 0 
# 95 percent confidence interval:
#  -1.803807   Inf 
# sample estimates:
# mean of x 
#17 

Hoping this helps!
Ted.


E-Mail: (Ted Harding) <[EMAIL PROTECTED]>
Fax-to-email: +44 (0)870 094 0861
Date: 09-Aug-08   Time: 22:28:34
-- XFMail --

__
R-help@r-project.org 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] Multivariate Regression with Weights

2008-08-09 Thread Arne Henningsen
On Saturday 09 August 2008 23:25, Arne Henningsen wrote:
> > -Original Message-
> > From: [EMAIL PROTECTED] [mailto:[EMAIL PROTECTED]
> > Sent: Monday, August 04, 2008 5:15 PM
> > To: Zhang Yanwei - Princeton-MRAm
> > Subject: RE: [R] Multivariate Regression with Weights
> >
> >   the systemfit package can do that and the documentation is quite nice
> > also.
>
> The systemfit() function of the systemfit package can use different weights
> for the different *equations* but NOT different weights for different
> *observations* (if argument "method" is "SUR", "WLS", "3SLS" or "W2SLS").
> I guess that you can divide the first equation by x_1 and the second
> equation by x_2 to obtain
>y_1^* = 1 + x_2^* + e_1^*
>y_2^* = 1 + x_1^* + e_2^*

Of course, I meant:
   y_1^* = b_11 + b_12  x_2^* + e_1^*
   y_2^* = b_22 + b_21  x_1^* + e_2^*

(Derived from
   y_1 = b_11  x_1 + b_12  x_2 + e_1
   y_2 = b_21  x_1 + b_22  x_2 + e_2
)

Sorry for my sloppiness.

Arne


> with
>y_1^* = y_1 / x_1
>x_2^* = x_2 / x_1
>e_1^* = e_1 / x_1
>y_2^* = y_2 / x_2
>x_1^* = x_1 / x_2
>e_2^* = e_2 / x_2
>
> In this case, you should have homoscedastic error terms. Hence, you can use
> systemfit to estimate the system of the modified equations.
>
> > On Mon, Aug 4, 2008 at  4:39 PM, Zhang Yanwei - Princeton-MRAm wrote:
> > > Hi all,
> > >   I'd like to fit a multivariate regression with the variance of the
> > > error term porportional to the predictors, like the WLS in the
> > > univariate case.
> > >   y_1~x_1+x_2
> > >   y_2~x_1+x_2
> > >   var(y_1)=x_1*sigma_1^2
>
> I guess that you mean var(e_1). Right?
>
> > >   var(y_2)=x_2*sigma_2^2
> > >   cov(y_1,y_2)=sqrt(x_1*x_2)*sigma_12^2
> > >
> > >  How can I specify this in R? Is there a corresponding function to the
> > > univariate specification lm(y~x,weights=x)?? Thanks.
>
> Unfortunately, systemfit has not (yet) a "weights" argument. (If you add
> this feature, I would be happy if you send me the code that I can include
> this feature in the official version.)
>
> Best wishes,
> Arne

-- 
Arne Henningsen
[EMAIL PROTECTED]
http://www.arne-henningsen.name


pgpFwCAx4Ot2M.pgp
Description: PGP signature
__
R-help@r-project.org 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 interpret t.test output

2008-08-09 Thread Felipe Carrillo
Hi Ted:
Thanks for your prompt reply and explanation. 
That's what I was wondering, why would one need to test mu=0 ,which is the 
t.test default. But reading Peter Dalgaard's book and looking at some examples 
online, I saw t.test being used like that; t.test(datasetname) with no other 
arguments.


> > t.test(fishlength)
> >  One Sample t-test
> > 
> > data:  fishlength
> > t = 30.1741, df = 13, p-value = 2.017e-13
> > alternative hypothesis: true mean is not equal to 0
> > 95 percent confidence interval:
> >  36.14141 41.71573
> > sample estimates:
> > mean of x
> >  38.92857
> > 
> > Thanks in advance for your help.
> 
> In terms of interpreting a statistical test, using your
> data,
> of the hypothesis that the mean length in the population is
> 0,
> the P-value of 0.2017 is very strong evidence
> indeed
> that the mean is not 0.
> 
> However, I do not know why you are asking the question. No
> test
> is needed. The length of any living fish, even while it is
> still
> in the egg, is greater than 0; and whatever population you
> have
> taken your sample from will have a mean length which is
> greater
> than 0.
> 
> That is not to say that the result of a t-test on any
> sample
> will necessarily give a significant result. You could have
> a
> small catch with lengths, say,
> 
>fishlengths <- c(2,4,9,20,50)
>t.test(fishlengths,mu=0)
> 
> # One Sample t-test
> # data:  fishlengths 
> # t = 1.9273, df = 4, p-value = 0.1262
> # alternative hypothesis: true mean is not equal to 0 
> # 95 percent confidence interval:
> #  -7.489442 41.489442 
> # sample estimates:
> # mean of x 
> #17 
> 
> And all you can conlude from that is that the sample, *in
> itself*,
> does not carry sufficient information to confirm what you
> know
> is true (i.e. mu > 0). Even the one-sided test of mu=0
> with alternative
> alt="greater" does not give a result significant
> at 5%:
> 
>  t.test(fishlengths,mu=0,alt="greater")
> 
> # One Sample t-test
> # data:  fishlengths 
> # t = 1.9273, df = 4, p-value = 0.0631
> # alternative hypothesis: true mean is greater than 0 
> # 95 percent confidence interval:
> #  -1.803807   Inf 
> # sample estimates:
> # mean of x 
> #17 
> 
> Hoping this helps!
> Ted.

__
R-help@r-project.org 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 interpret t.test output

2008-08-09 Thread Felipe Carrillo
I take it back, Peter Dalgaard's book uses t.test with mu=7725 and no mu=0. I 
got the script online.

Hi Ted:
Thanks for your prompt reply and explanation. 
That's what I was wondering, why would one need to test mu=0 ,which is the 
t.test default. But reading Peter Dalgaard's book and looking at some examples 
online, I saw t.test being used like that; t.test(datasetname) with no other 
arguments.


> > t.test(fishlength)
> >  One Sample t-test
> > 
> > data:  fishlength
> > t = 30.1741, df = 13, p-value = 2.017e-13
> > alternative hypothesis: true mean is not equal to 0
> > 95 percent confidence interval:
> >  36.14141 41.71573
> > sample estimates:
> > mean of x
> >  38.92857
> > 
> > Thanks in advance for your help.
> 
> In terms of interpreting a statistical test, using your
> data,
> of the hypothesis that the mean length in the population is
> 0,
> the P-value of 0.2017 is very strong evidence
> indeed
> that the mean is not 0.
> 
> However, I do not know why you are asking the question. No
> test
> is needed. The length of any living fish, even while it is
> still
> in the egg, is greater than 0; and whatever population you
> have
> taken your sample from will have a mean length which is
> greater
> than 0.
> 
> That is not to say that the result of a t-test on any
> sample
> will necessarily give a significant result. You could have
> a
> small catch with lengths, say,
> 
>fishlengths <- c(2,4,9,20,50)
>t.test(fishlengths,mu=0)
> 
> # One Sample t-test
> # data:  fishlengths 
> # t = 1.9273, df = 4, p-value = 0.1262
> # alternative hypothesis: true mean is not equal to 0 
> # 95 percent confidence interval:
> #  -7.489442 41.489442 
> # sample estimates:
> # mean of x 
> #17 
> 
> And all you can conlude from that is that the sample, *in
> itself*,
> does not carry sufficient information to confirm what you
> know
> is true (i.e. mu > 0). Even the one-sided test of mu=0
> with alternative
> alt="greater" does not give a result significant
> at 5%:
> 
>  t.test(fishlengths,mu=0,alt="greater")
> 
> # One Sample t-test
> # data:  fishlengths 
> # t = 1.9273, df = 4, p-value = 0.0631
> # alternative hypothesis: true mean is greater than 0 
> # 95 percent confidence interval:
> #  -1.803807   Inf 
> # sample estimates:
> # mean of x 
> #17 
> 
> Hoping this helps!
> Ted.

__
R-help@r-project.org 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] Converting nested "for" loops to an "apply" function(s)

2008-08-09 Thread Kurt Newman

Hello,

I would like to know more about how to use the "apply" family and have 
attempted to convert nested "for" loops in example code from Contributed 
Documentation ("The Friendly Beginners' R Course” by Toby Marthews (ZIP, 
2007-03-01)") to an "apply" function(s).  The relevant code is:

distances=c(51,65,175,196,197,125,10,56)#distances of 8 houses from the 
town centre in m
bearings=c(10,8,210,25,74,128,235,335)  #bearings of the houses in degrees

xpos=distances*sin(bearings*pi/180) #in sin and cos the argument 
MUST be in radians
ypos=distances*cos(bearings*pi/180) 

numpoints=length(distances)
nnd=rep(sqrt(2*400*400),times=numpoints)#start with the maximum 
possible distance
for (i in 1:numpoints) {
 for (j in 1:numpoints) {
  if (i!=j) {   
   diffx=abs(xpos[i]-xpos[j])
   diffy=abs(ypos[i]-ypos[j])
   nd=sqrt((diffx^2)+(diffy^2))
   if (nd
__
R-help@r-project.org 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 interpret t.test output

2008-08-09 Thread Peter Dalgaard

Felipe Carrillo wrote:

Hi Ted:
Thanks for your prompt reply and explanation. 
That's what I was wondering, why would one need to test mu=0 ,which is the t.test default. But reading Peter Dalgaard's book and looking at some examples online, I saw t.test being used like that; t.test(datasetname) with no other arguments

Umm, in ISwR, I  can see

t.test(bmi, mu=22.5)
t.test(bmi, mu=22.5)$p.value
t.test(daily.intake,mu=7725)
t.test(daily.intake,mu=7725)
t.test(expend~stature)
t.test(expend~stature, var.equal=T)
t.test(pre, post, paired=T)
t.test(pre, post) #WRONG!
t.test(log10(diameter)~glucose)

None of those are of the form t.test(mydata). 

You might occasionally want to  test for mu=0, for instance in 
t.test(post-pre), which is  just another way of doing the paired t test. 
mu=0 is the default because zero is the only value that "sticks out" as 
a potentially  hypothesized value (in particular, it is stable to 
scaling of the data).


If you see cases where mu is not specified and the default 0 is not 
interesting as a null hypothesis, then maybe the author was not 
interested the test at all (e.g. the confidence interval is still useful).


--
  O__   Peter Dalgaard Øster Farimagsgade 5, Entr.B
 c/ /'_ --- Dept. of Biostatistics PO Box 2099, 1014 Cph. K
(*) \(*) -- University of Copenhagen   Denmark  Ph:  (+45) 35327918
~~ - ([EMAIL PROTECTED])  FAX: (+45) 35327907

__
R-help@r-project.org 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] effective matrix subset

2008-08-09 Thread jgarcia
It seems that this solution provided by Dan (and also available in
SPoetry; I'm sorry I didn't notice it) is the fastest and simplest. I was
using a more standard approach:

V <- t(A)[(0:(nrow(A)-1))*ncol(A)+X],

That wasn't bad, but I was confident that you, R gurus, could outperform
this. This is clearly a much better solution.

Thanks all, and best wishes,
Javier
--




> on 08/09/2008 06:52 AM Dan Davison wrote:
>> On Sat, Aug 09, 2008 at 06:29:59AM -0500, Marc Schwartz wrote:
>>> on 08/09/2008 06:01 AM [EMAIL PROTECTED] wrote:
 Hi;
 If we have a matrix A, and a vector X, where length(X)=nrow(A), and X
 contains a wanted column for each row in A, in row ascending order.
 How
 would be the most effective way to extract the desired vector V (with
 length(V)=nrow(A))?
>>>
>>> A <- matrix(1:20, 4, 5)
>>>
 A
>>>  [,1] [,2] [,3] [,4] [,5]
>>> [1,]159   13   17
>>> [2,]26   10   14   18
>>> [3,]37   11   15   19
>>> [4,]48   12   16   20
>>>
>>>
>>> # Create an arbitrary set of indices, one for each row in A
>>> X <- c(2, 5, 1, 4)
>>>
 X
>>> [1] 2 5 1 4
>>>
>>>
>>> Presumably you want:
>>>
>>> V <- c(A[1, 2], A[2, 5], A[3, 1], A[4, 4])
>>>
 V
>>> [1]  5 18  3 16
>>>
>>>
>>> If so, then:
>>>
 sapply(seq(nrow(A)), function(i) A[i, X[i]])
>>> [1]  5 18  3 16
>>
>> Or
>>
>>> A[cbind(seq(nrow(A)), X)]
>> [1]  5 18  3 16
>>
>> Dan
>
> Better (and faster) solution Dan.
>
> I can't blame the lack of coffee on missing that one this morning. I
> have had a full pot already over the past 6 hours, working on shifting
> my internal clock and getting ready to begin my journey to Dortmund
> later tonight...
>
> Safe travels to all who are going.
>
> Marc
>
>

__
R-help@r-project.org 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 repress the annoying complains from X window system

2008-08-09 Thread dusa.adrian

Dear Marc,


Marc Schwartz wrote:
> 
> On Fri, 2007-05-18 at 11:25 -0400, Hao Liu wrote:
>> Dear All:
>> 
>> I am running some GUI functions in linux environment, they runs fine, 
>> however I constantly get this kind of message in R console:
>> 
>> Warning: X11 protocol error: BadWindow (invalid Window parameter)
>> 
>> Is there a way to repress it? Or am I doing something wrong here.. it 
>> does not interfere with the running of fucntion though.
>> 
>> Thanks
>> Hao
> 
> Upgrade your version of R.
> 
> You have not provided sufficient details, but if I had to guess, you are
> either running RCmdr or using other tcl/tk based widgets.
> 
> If correct, the error message that you are seeing was fixed back in R
> 2.4.0:
> 
> o The X11() device no longer produces (apparently spurious)
>   'BadWindow (invalid Window parameter)' warnings when run from
>   Rcmdr.
> 

I'm digging on a past message, where you said that message was fixed in R
2.4.0.
However, I'm using R 2.7.1 (under Kubuntu Hardy Heron, with tcl and tk
version 8.4) and still experience the same message with Rcmdr.
For reproduction, this happens when using (before loading Rcmdr):
options(Rcmdr=list(console.output=TRUE)).

In addition to this error message, I also get sometimes this one (repeated
17 times at once):

Warning in structure(.External("dotTclObjv", objv, PACKAGE = "tcltk"), class
= "tclObj") :
  X11 protocol error: BadWindow (invalid Window parameter)

> R.version
   _
platform   i486-pc-linux-gnu
arch   i486
os linux-gnu
system i486, linux-gnu
status
major  2
minor  7.1
year   2008
month  06
day23
svn rev45970
language   R
version.string R version 2.7.1 (2008-06-23)

Thanks in advance for any hint,
Adrian

-- 
View this message in context: 
http://www.nabble.com/How-to-repress-the-annoying-complains-from-X-window-system-tp10684882p18909311.html
Sent from the R help mailing list archive at Nabble.com.

__
R-help@r-project.org 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] box-plot without a box surrounding the plot

2008-08-09 Thread Jörg Groß

Hi,

I want to draw a boxplot and I want to get rid of the surrounding box.
So that there are only the axis lines left on the bottom and on the  
left.


I tried a litte bit with box() but I dont get it the way I want.

Can somebody help me out?

__
R-help@r-project.org 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] box-plot without a box surrounding the plot

2008-08-09 Thread stephen sefick
par(bty="l")
boxplot(count ~ spray, data = InsectSprays, col = "lightgray")
#is that what you want?

Stephen
On Sat, Aug 9, 2008 at 8:12 PM, Jörg Groß <[EMAIL PROTECTED]> wrote:
> Hi,
>
> I want to draw a boxplot and I want to get rid of the surrounding box.
> So that there are only the axis lines left on the bottom and on the left.
>
> I tried a litte bit with box() but I dont get it the way I want.
>
> Can somebody help me out?
>
> __
> R-help@r-project.org 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.
>



-- 
Let's not spend our time and resources thinking about things that are
so little or so large that all they really do for us is puff us up and
make us feel like gods. We are mammals, and have not exhausted the
annoying little problems of being mammals.

-K. Mullis

__
R-help@r-project.org 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] possible problem with rgdal

2008-08-09 Thread Naiara Pinto
Hello all,

I recently installed rgdal 0.5.24-1 (kyngchaos framework) and I am having
trouble making sense of the row, col information provided.

> a = new("GDALDataset", "dummy.tif")
> GDALinfo("dummy.tif")
rows420
columns 660
bands   1
ll.x-55.5
ll.y-14.5
res.x   0.00834
res.y   0.00834
oblique.x   0
oblique.y   0
driver  GTiff
projection  +proj=longlat +ellps=WGS84 +datum=WGS84 +no_defs
filedummy.tif
> nrow(a)
[1] 420
> ncol(a)
[1] 660
> dummy = getRasterData(a)
> nrow(dummy)
[1] 660
> ncol(dummy)
[1] 420
> test = matrix(0, nrow=420, ncol=660)
> putRasterData(a,test)
Error in putRasterData(a, test) :
GDAL Error 5: Access window out of range in RasterIO().  Requested
(0,0) of size 420x660 on raster of 660x420.

Also, it appears that the lower left coordinates given by calling GDALinfo
are actually the upper left coordinates.

from command line GDAL:
sardinha$ gdalinfo dummy.tif
Driver: GTiff/GeoTIFF
Files: dummy.tif
Size is 660, 420
Coordinate System is:
GEOGCS["WGS 84",
DATUM["WGS_1984",
SPHEROID["WGS 84",6378137,298.2572235630016,
AUTHORITY["EPSG","7030"]],
AUTHORITY["EPSG","6326"]],
PRIMEM["Greenwich",0],
UNIT["degree",0.0174532925199433],
AUTHORITY["EPSG","4326"]]
Origin = (-55.500,-14.500)
Pixel Size = (0.00833767951,-0.00833767951)
Metadata:
  AREA_OR_POINT=Area
Image Structure Metadata:
  INTERLEAVE=BAND
Corner Coordinates:
Upper Left  ( -55.500, -14.500) ( 55d30'0.00"W, 14d30'0.00"S)
Lower Left  ( -55.500, -18.002) ( 55d30'0.00"W, 18d 0'0.00"S)
Upper Right ( -49.997, -14.500) ( 50d 0'0.00"W, 14d30'0.00"S)
Lower Right ( -49.997, -18.002) ( 50d 0'0.00"W, 18d 0'0.00"S)
Center  ( -52.749, -16.251) ( 52d45'0.00"W, 16d15'0.00"S)

I am not sure this is a problem, but I don't remember seeing this before.
Does anybody know if there is an older version that gives the correct
information?

Thanks,

Naiara.

[[alternative HTML version deleted]]

__
R-help@r-project.org 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] box-plot without a box surrounding the plot

2008-08-09 Thread Jörg Groß




par(bty="l")
boxplot(count ~ spray, data = InsectSprays, col = "lightgray")
#is that what you want?

Stephen



Thanks! That's what I needed.




On Sat, Aug 9, 2008 at 8:12 PM, Jörg Groß <[EMAIL PROTECTED]>  
wrote:

Hi,

I want to draw a boxplot and I want to get rid of the surrounding  
box.
So that there are only the axis lines left on the bottom and on the  
left.


I tried a litte bit with box() but I dont get it the way I want.

Can somebody help me out?

__
R-help@r-project.org 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.





--
Let's not spend our time and resources thinking about things that are
so little or so large that all they really do for us is puff us up and
make us feel like gods. We are mammals, and have not exhausted the
annoying little problems of being mammals.

-K. Mullis


__
R-help@r-project.org 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] Converting nested "for" loops to an "apply" function(s)

2008-08-09 Thread Kurt Newman

Resending.  Previous message was truncated.  Sorry for possible confusion.


> From: [EMAIL PROTECTED]
> To: r-help@r-project.org
> Date: Sat, 9 Aug 2008 18:25:47 -0400
> Subject: [R] Converting nested "for" loops to an "apply" function(s)
> 
> 
> Hello,
> 
> I would like to know more about how to use the "apply" family and have 
> attempted to convert nested "for" loops in example code from Contributed 
> Documentation ("The Friendly Beginners' R Course” by Toby Marthews (ZIP, 
> 2007-03-01)") to an "apply" function(s).  The relevant code is:
> 
> distances=c(51,65,175,196,197,125,10,56)  #distances of 8 houses from the 
> town centre in m
> bearings=c(10,8,210,25,74,128,235,335)#bearings of the houses in 
> degrees
> 
> xpos=distances*sin(bearings*pi/180)   #in sin and cos the argument 
> MUST be in radians
> ypos=distances*cos(bearings*pi/180)   
> 
> numpoints=length(distances)
> nnd=rep(sqrt(2*400*400),times=numpoints)  #start with the maximum 
> possible distance
> for (i in 1:numpoints) {  
>  for (j in 1:numpoints) {
>   if (i!=j) { 
>diffx=abs(xpos[i]-xpos[j])
>diffy=abs(ypos[i]-ypos[j])
>nd=sqrt((diffx^2)+(diffy^2))
>if (nd < nnd[i]) {nnd[i]=nd}
  }
 }
}
print(data.frame(xpos,ypos,nnd))

My attempts to convert the nested "for" loops to an "apply" function(s) have 
not been successful.  I would like to know how to convert the code to increase 
my knowledge of R programming and to evaluate operational efficiency of the 
different strategies.

Thank you in advance for your comments / suggestions.

Kurt Newman



> __
> R-help@r-project.org 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@r-project.org 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] Scripting - query

2008-08-09 Thread Gareth Campbell
I have a vector:
alleles.present<-c("D3", "D16", ... )

The alleles present changes given the case I'm dealing with - i.e. either
all of the alleles I use for my calculations are present, or some of them.

Depending on what alleles are present, I need to make matrices and do
calculations on those alleles present and completely disregard any formula
or other use of the alleles not present.

I'm trying to figure out the best way to do this.

Basically I'm trying to do if() commands (with no success so far) to allow
me to query the alleles.present for the presence of each allele I use and
then let dictate which formula to use etc...

Does anyone have a good way to do this?  I've been fiddling with grep()
etc... but I can't get it to do what I need!!  Very frustrating.

Thanks very much

-- 
Gareth Campbell
PhD Candidate
The University of Auckland

P +649 815 3670
M +6421 256 3511
E [EMAIL PROTECTED]
[EMAIL PROTECTED]

[[alternative HTML version deleted]]

__
R-help@r-project.org 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] Reading large datasets and fitting logistic models in R

2008-08-09 Thread Pradheep K E
Hi R-experts,

Does anyone have experience using R for handling large scale data (millions
of rows, hundreds or thousands of features)?

What is the largest size of data that anyone has used with glm?

Also, is there a library to read data in sparse data format (like SVMlight
format)?

Thanks
Pradheep

[[alternative HTML version deleted]]

__
R-help@r-project.org 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] instal tar.gz package on windows

2008-08-09 Thread Andrew Harris
For cross-building a mature package, I stumbled across a convenient 
utility provided by the R Windows maintainer:

_http://win-builder.r-project.org/

_This takes a little while, so for iteration it is much better to 
install the Rtools, but for a package or two it is great.


Best,
Andrew Harris, Univ. Maryland

__
R-help@r-project.org 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] mean-plot (add residuals)

2008-08-09 Thread Jörg Groß

Hi,

I want to plot the mean of a variable and add the standard deviation  
as a line going above and below the mean (like the whiskers in a  
boxplot) but I don't know how to add these residual-lines.


Is there an easy way to do that or do I have to use the line()- 
function to create them on my own?



Thanks for help!

__
R-help@r-project.org 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] instal tar.gz package on windows

2008-08-09 Thread Prof Brian Ripley

In apparently a late response to

https://stat.ethz.ch/pipermail/r-help/2008-July/169578.html

on Sat, 9 Aug 2008, Andrew Harris wrote:

For cross-building a mature package, I stumbled across a convenient utility 
provided by the R Windows maintainer:

_http://win-builder.r-project.org/


It is in both the 'Writing R Extensions' and 'R Installation and 
Administration' manuals, so do we take it you have not yet 'stumbled 
across' those?


There's lots of useful information in the R manuals!

_This takes a little while, so for iteration it is much better to install the 
Rtools, but for a package or two it is great.


Best,
Andrew Harris, Univ. Maryland

__
R-help@r-project.org 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.


--
Brian D. Ripley,  [EMAIL PROTECTED]
Professor of Applied Statistics,  http://www.stats.ox.ac.uk/~ripley/
University of Oxford, Tel:  +44 1865 272861 (self)
1 South Parks Road, +44 1865 272866 (PA)
Oxford OX1 3TG, UKFax:  +44 1865 272595

__
R-help@r-project.org 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] Reading large datasets and fitting logistic models in R

2008-08-09 Thread Prof Brian Ripley

See also bigglm() in package biglm.

On Sat, 9 Aug 2008, Pradheep K E wrote:


Hi R-experts,

Does anyone have experience using R for handling large scale data (millions
of rows, hundreds or thousands of features)?

What is the largest size of data that anyone has used with glm?


I've used 700,000 rows and about 100 cols, but it was 4 years ago and we 
have more memory now.  It matters if the 'features' are numeric or 
categorical, as the latter can expand to many columns in the model matrix.


As a rough guide, expect to need 200x as much memory in bytes as nrows x 
ncols.  Using glm.fit will be more efficient (I've just tested 100,000 x 
100 which used 1.2Gb).



Also, is there a library to read data in sparse data format (like SVMlight
format)?


You mean *store* data in a sparse format when read in?  I'm not sure of 
the relevance, but look at the function method for bigglm for a way to 
avoid even doing that. If the data are numeric there are at least three 
sparse-matrix packages on CRAN.


Ultimately R's code such as glm() is designed for flexibility and to do 
interesting things with the fit: for really large problems you will do 
better to write a specialized fitting routine.  bigglm() is an

intermediate position.

There's also the question of whether there are any interesting homogeneous 
datasets of this sort of size.  Often doing analyses on subsets and a 
meta-analysis is a much more insightful approach (as it was in our 
problem: we split on one of the categorical explanatory variables).



Thanks
Pradheep

[[alternative HTML version deleted]]

__
R-help@r-project.org 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.



--
Brian D. Ripley,  [EMAIL PROTECTED]
Professor of Applied Statistics,  http://www.stats.ox.ac.uk/~ripley/
University of Oxford, Tel:  +44 1865 272861 (self)
1 South Parks Road, +44 1865 272866 (PA)
Oxford OX1 3TG, UKFax:  +44 1865 272595

__
R-help@r-project.org 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] One way ANOVA

2008-08-09 Thread Angelo Scozzarella

Hi,

There is  a command or a function for One way ANOVA if I know the  
numbers of subjects for each group, the mean and di standard deviation  
for each group?


Thanks


Angelo Scozzarella

__
R-help@r-project.org 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.