Re: [R] soil.texture() function with bug?
Hi Patrick, Not sure what the problem is from your email. Can you send a code example which reproduces the error? Make sure you mention the version of R you are using! Also send your question/reply to/cc R-help as well. Then the rest of the world is there to help you too. Greetings, Sander. Patrick Kuss wrote: Dear Sander Oom and Jim Lemon, thanks for putting the soils.texture() function into R. However, for whatever reason I am not able to display the triangle correctly. Each of the 27 tick labels shows as c(10,20,30,40,50,60,70,80,90) and thus basically cover the whole triangle. Plotting points and adding graphical paramters like 'col.symbols' works fine. I am also able to plot my soil data using triax.plot() without any problems, but of course, the nice feature of soil type areas is not available. Do you have any insights? Thanks a lot and cheers from Alaska Patrick __ 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] Data frame referencing?
Dear R users, When you do: x - rnorm(10) y - rnorm(10) z - rnorm(10) a - data.frame(x,y,z) a$x [1] 1.37821893 0.21152756 -0.55453182 -2.10426048 -0.08967880 0.03712110 [7] -0.80592149 0.07413450 0.15557671 1.22165341 Why does this not work: a[a$y0.5,y] -1 Error in [-.data.frame(`*tmp*`, a$y 0.5, y, value = 1) : only 0's may be mixed with negative subscripts While this works: a[a$y0.5,2] -1 a x y z 1 1.37821893 -1.0887363 1.7340522 2 0.21152756 -0.7256467 -1.3165373 3 -0.55453182 1.000 -2.1116072 4 -2.10426048 -0.4898596 -1.5863823 5 -0.08967880 1.000 -0.9139706 6 0.03712110 1.000 -1.3004970 7 -0.80592149 -0.7004193 -0.1958059 8 0.07413450 1.000 -1.3574303 9 0.15557671 -0.3335407 -2.1991236 10 1.22165341 1.000 -0.7576708 For a complex loop I would prefer to reference the right colomn by name, not by number! Now, when the colomns change, I need to check my code to make sure that the right colomns are referenced. Suggestions much appreciated! Thanks in advance, Sander. __ 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] Data frame referencing?
Dear Gabor and Dimitris, Simple, once you know! It is these little exceptions on the R notation that get me stuck. Now I am on the loose again! Thanks, Sander. Dimitris Rizopoulos wrote: you need to use quotes, i.e., a[a$y 0.5, y] - 1 you can also use a$y[a$y 0.5] - 1 I hope it helps. Best, Dimitris Dimitris Rizopoulos Ph.D. Student Biostatistical Centre School of Public Health Catholic University of Leuven Address: Kapucijnenvoer 35, Leuven, Belgium Tel: +32/(0)16/336899 Fax: +32/(0)16/337015 Web: http://med.kuleuven.be/biostat/ http://www.student.kuleuven.be/~m0390867/dimitris.htm - Original Message - From: Sander Oom [EMAIL PROTECTED] To: r-help@stat.math.ethz.ch Sent: Friday, August 04, 2006 1:48 PM Subject: [R] Data frame referencing? Dear R users, When you do: x - rnorm(10) y - rnorm(10) z - rnorm(10) a - data.frame(x,y,z) a$x [1] 1.37821893 0.21152756 -0.55453182 -2.10426048 -0.08967880 0.03712110 [7] -0.80592149 0.07413450 0.15557671 1.22165341 Why does this not work: a[a$y0.5,y] -1 Error in [-.data.frame(`*tmp*`, a$y 0.5, y, value = 1) : only 0's may be mixed with negative subscripts While this works: a[a$y0.5,2] -1 a x y z 1 1.37821893 -1.0887363 1.7340522 2 0.21152756 -0.7256467 -1.3165373 3 -0.55453182 1.000 -2.1116072 4 -2.10426048 -0.4898596 -1.5863823 5 -0.08967880 1.000 -0.9139706 6 0.03712110 1.000 -1.3004970 7 -0.80592149 -0.7004193 -0.1958059 8 0.07413450 1.000 -1.3574303 9 0.15557671 -0.3335407 -2.1991236 10 1.22165341 1.000 -0.7576708 For a complex loop I would prefer to reference the right colomn by name, not by number! Now, when the colomns change, I need to check my code to make sure that the right colomns are referenced. Suggestions much appreciated! Thanks in advance, Sander. __ 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. Disclaimer: http://www.kuleuven.be/cwis/email_disclaimer.htm __ 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] Conditional operation on columns in data frame
Dear R-users,
I need to do an SQL like, conditional, operation on a data frame using
an ifelse construction. However I can not get it to work.
Example:
x1 - rnorm(10)
x2 - rnorm(10)
x3 - rnorm(10)
x3 - NA
y - cbind(x1,x2,x3)
y
x1 x2 x3
[1,] -0.56927780 -0.30952274 NA
[2,] 0.16355087 0.05911772 NA
[3,] -0.21899354 2.04583752 NA
[4,] 0.06368076 1.11661608 NA
[5,] -1.30249878 -0.63354373 NA
[6,] 0.04842365 1.47591928 NA
[7,] -0.32364275 -0.62201121 NA
[8,] -0.70427823 -0.15485223 NA
[9,] -0.39563916 2.23504977 NA
[10,] -1.24102239 -0.40991140 NA
Now I want to fill a new column with values from either x2 or x3,
depending on the value in x1. I thought of something like this:
y$x4 - ifelse(y$x10,x2,x3)
y$x4 - apply(y,1, function(x) {ifelse(x$x10,x$x2,x$x3) })
But obviously this did not give the right result.
Any suggestions?
Thanks in advance,
Sander.
__
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
Re: [R] Conditional operation on columns in data frame
Hi Hong Ooi,
Thanks for your reply! I tried this option as well, but the result is
not as I expect:
y$x4 - ifelse(y$x10,y$x2,y$x3)
y
$x4
logical(0)
class(y)
[1] list
Any more ideas?
Thanks,
Sander.
Hong Ooi wrote:
___
Note: This e-mail is subject to the disclaimer contained at the bottom of
this message.
___
x2 and x3 are not separate objects, but columns of y. So you need to specify
that.
y$x4 - ifelse(y$x1 0, y$x2, y$x3)
From: [EMAIL PROTECTED] on behalf of Sander Oom
Sent: Sat 11/03/2006 11:23 PM
To: r-help@stat.math.ethz.ch
Subject: [R] Conditional operation on columns in data frame
Dear R-users,
I need to do an SQL like, conditional, operation on a data frame using
an ifelse construction. However I can not get it to work.
Example:
x1 - rnorm(10)
x2 - rnorm(10)
x3 - rnorm(10)
x3 - NA
y - cbind(x1,x2,x3)
y
x1 x2 x3
[1,] -0.56927780 -0.30952274 NA
[2,] 0.16355087 0.05911772 NA
[3,] -0.21899354 2.04583752 NA
[4,] 0.06368076 1.11661608 NA
[5,] -1.30249878 -0.63354373 NA
[6,] 0.04842365 1.47591928 NA
[7,] -0.32364275 -0.62201121 NA
[8,] -0.70427823 -0.15485223 NA
[9,] -0.39563916 2.23504977 NA
[10,] -1.24102239 -0.40991140 NA
Now I want to fill a new column with values from either x2 or x3,
depending on the value in x1. I thought of something like this:
y$x4 - ifelse(y$x10,x2,x3)
y$x4 - apply(y,1, function(x) {ifelse(x$x10,x$x2,x$x3) })
But obviously this did not give the right result.
Any suggestions?
Thanks in advance,
Sander.
__
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
___
The information transmitted in this message and its attach...{{dropped}}
__
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
Re: [R] How to repeat code snippet for several variables in a data frame?
Thanks for the very useful tips! Now I have enough round and square bracket and other tricks to wrap up the function! The double square bracket trick in test[[varname]] is golden! Thanks again, Sander. __ 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
[R] How to repeat code snippet for several variables in a data frame?
Dear all, I have a data frame containing the results of an experiment. Like this: a-seq(1,4,by=1) b-seq(1,2,by=1) test-expand.grid(b,a,a) colnames(test)-c(replicates,bins, groups) test$abc - rnorm(32) test$def - rnorm(32) test$ghi - rnorm(32) test The following code snippet aggregates the data for one variable and then draws a plot: tmp - aggregate(test$abc, list( test$bins, test$groups), mean) colnames(tmp) - c(bins, groups, abc) tmp #pltName - paste(line_datGrassChem_, abc, .eps, sep=) #postscript(pltName) labs - c(-15/-9,-9/-6,-6/-3,-3/0) sps - trellis.par.get(superpose.symbol) sps$pch - 1:4 trellis.par.set(superpose.symbol, sps) xyplot( abc ~ bins, data = tmp, groups = groups, type = list(p, l), scales = list(x = list(labels=labs)), layout = c(1,1), key = list(columns = 4, text = list(paste(unique(tmp$groups))), points = Rows(sps, 1:4) ) ) #dev.off() How can I wrap R code around this code snippet, such that I can repeat the same code snippet for all other variables in the data frame (i.e. def, ghi, etc.). Thanks for your suggestions! Sander. -- Dr Sander P. Oom Animal, Plant and Environmental Sciences, University of the Witwatersrand Private Bag 3, Wits 2050, South Africa Tel (work) +27 (0)11 717 64 04 Tel (home) +27 (0)18 297 44 51 Fax +27 (0)18 299 24 64 Email [EMAIL PROTECTED] Web www.oomvanlieshout.net/sander __ 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
Re: [R] some suggestion
Shame nobody has included this function in their package! Would be useful to have as a standard function! Sander. Dimitris Rizopoulos wrote: take a look at this function by Kevin Wright RSiteSearch(sort.data.frame) Best, Dimitris Dimitris Rizopoulos Ph.D. Student Biostatistical Centre School of Public Health Catholic University of Leuven Address: Kapucijnenvoer 35, Leuven, Belgium Tel: +32/16/336899 Fax: +32/16/337015 Web: http://www.med.kuleuven.ac.be/biostat/ http://www.student.kuleuven.ac.be/~m0390867/dimitris.htm - Original Message - From: ronggui [EMAIL PROTECTED] To: r-help@stat.math.ethz.ch Sent: Tuesday, June 14, 2005 4:28 PM Subject: [R] some suggestion it seems R has no function to sort the data.frame according to some variable(s),though we can do these by order() and indexing.but why not make sort() the a generic function,making it can sort vector and other type of objects? maybe this is a silly suggestion,but i think it is quite usefull. 2005-06-14 -- Deparment of Sociology Fudan University Blog:www.sociology.yculblog.com __ 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 __ 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 -- Dr Sander P. Oom Animal, Plant and Environmental Sciences, University of the Witwatersrand Private Bag 3, Wits 2050, South Africa Tel (work) +27 (0)11 717 64 04 Tel (home) +27 (0)18 297 44 51 Fax +27 (0)18 299 24 64 Email [EMAIL PROTECTED] Web www.oomvanlieshout.net/sander __ 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
Re: [R] Replacing for loop with tapply!?
Dear Adaikalavan,
Your solution (the second function) is definitely the most elegant and
generic solution of all replies in this discussion. Robust for missing
values and flexible to allow as many calculations as desired! It is so
clear, I even managed to hack it (of course also thanks to the new
insight from all the other posts)!
As the data consists of weather stations in rows and days in columns, I
have adapted the function to work on rows instead of columns. Did not
manage to get the results directly into the right rows/cols layout, so a
transpose (t) is still required. However this seems instant, so does not
mean a reduction in speed! Calculating proportions is now a snip!!
Thanks for you help,
Sander.
### simulate data
set.seed(1)# for reproducibility
mat - matrix(sample(-15:50, 15 * 10, TRUE), 15, 10)
mat[ mat 45 ] - NA # create some missing values
mat[ 9, ] - NA # station 9's data is completely missing
mat
find.stats - function( data, threshold ){
n - length(threshold)
excess - numeric( n )
out- matrix( ncol=nrow(data), nrow=(n + 2) ) # initialise
good - which( apply( data, 1, function(x) !all(is.na(x)) ) )
# rows that are not completely missing
out[ ,good ] - apply( data[ good, ], 1, function(x){
m - max( x, na.rm=T )
# determine maximum value per row
c - length(x[!is.na(x)])
# determine number of non-missing values
for(i in 1:n){ excess[i] - sum( x threshold[i], na.rm=TRUE
)/length(x[!is.na(x)]) }
# calc proportion of non-missing values over multiple thresholds
return( c(m, c, excess) )
} )
rownames(out) - c( TmpMax, Count, paste(Over, threshold, sep=) )
colnames(out) - rownames(data) # name of the stations
return( t(out) )
}
lstTemps=c(37,39,41,43)
tmp - find.stats( mat, lstTemps )
tmp
Adaikalavan Ramasamy wrote:
OK, so you want to find some summary statistics for each column, where
some columns could be completely missing.
Writing a small wrapper should help. When you use apply(), you are
actually applying a function to every column (or row). First, let us
simulate a dataset with 15 days/rows and 10 stations/columns
### simulate data
set.seed(1)# for reproducibility
mat - matrix(sample(-15:50, 15 * 10, TRUE), 15, 10)
mat[ mat 45 ] - NA # create some missing values
mat[ ,9 ] - NA # station 9's data is completely missing
Here are two example of such wrappers :
find.stats1 - function( data, threshold=c(37,39,41) ){
n - length(threshold)
out - matrix( nrow=(n + 1), ncol=ncol(data) ) # initialise
out[1, ] - apply(data, 2, function(x)
ifelse( all(is.na(x)), NA, max(x, na.rm=T) ))
for(i in 1:n) out[ i+1, ] - colSums( data threshold[i], na.rm=T )
rownames(out) - c( daily_max, paste(above, threshold, sep=_) )
colnames(out) - rownames(data) # name of the stations
return( out )
}
find.stats2 - function( data, threshold=c(37,39,41) ){
n - length(threshold)
excess - numeric( n )
out- matrix( nrow=(n + 1), ncol=ncol(data) ) # initialise
good - which( apply( data, 2, function(x) !all(is.na(x)) ) )
# colums that are not completely missing
out[ , good] - apply( data[ , good], 2, function(x){
m - max( x, na.rm=T )
for(i in 1:n){ excess[i] - sum( x threshold[i], na.rm=TRUE ) }
return( c(m, excess) )
} )
rownames(out) - c( daily_max, paste(above, threshold, sep=_) )
colnames(out) - rownames(data) # name of the stations
return( out )
}
find.stats1( mat )
[,1] [,2] [,3] [,4] [,5] [,6] [,7] [,8] [,9] [,10]
daily_max 44 42 39 41 45 43 42 45 NA42
above_37 212132210 1
above_39 210132110 1
above_41 210022110 1
find.stats2( mat )
[,1] [,2] [,3] [,4] [,5] [,6] [,7] [,8] [,9] [,10]
daily_max 44 42 39 41 45 43 42 45 NA42
above_37 21213221 NA 1
above_39 21013211 NA 1
above_41 21002211 NA 1
On my laptop 'find.stats1' and 'find.stats2' (which is more flexible)
takes 7 and 6 seconds respectively to execute on a dataset with 1
stations and 365 days.
Regards, Adai
On Fri, 2005-06-10 at 20:05 +0200, Sander Oom wrote:
Dear all,
Dimitris and Andy, thanks for your great help. I have progressed to the
following code which runs very fast and effective:
mat - matrix(sample(-15:50, 15 * 10, TRUE), 15, 10)
mat[mat45] - NA
mat-NA
mat
temps - c(35, 37, 39)
ind - rbind(
t(sapply(temps, function(temp)
rowSums(mat temp, na.rm=TRUE) )),
rowSums(!is.na(mat), na.rm=FALSE),
apply(mat, 1, max, na.rm=TRUE))
ind - t(ind)
ind
However, some weather stations have missing values for the whole
[R] Replacing for loop with tapply!?
Dear all,
We have a large data set with temperature data for weather stations
across the globe (15000 stations).
For each station, we need to calculate the number of days a certain
temperature is exceeded.
So far we used the following S code, where mat88 is a matrix containing
rows of 365 daily temperatures for each of 15000 weather stations:
m - 37
n - 2
outmat88 - matrix(0, ncol = 4, nrow = nrow(mat88))
for(i in 1:nrow(mat88)) {
# i - 3
row1 - as.data.frame(df88[i, ])
temprow37 - select.rows(row1, row1 m)
temprow39 - select.rows(row1, row1 m + n)
temprow41 - select.rows(row1, row1 m + 2 * n)
outmat88[i, 1] - max(row1, na.rm = T)
outmat88[i, 2] - count.rows(temprow37)
outmat88[i, 3] - count.rows(temprow39)
outmat88[i, 4] - count.rows(temprow41)
}
outmat88
We have transferred the data to a more potent Linux box running R, but
still hope to speed up the code.
I know a for loop should be avoided when looking for speed. I also know
the answer is in something like tapply, but my understanding of these
commands is still to limited to see the solution. Could someone show me
the way!?
Thanks in advance,
Sander.
--
Dr Sander P. Oom
Animal, Plant and Environmental Sciences,
University of the Witwatersrand
Private Bag 3, Wits 2050, South Africa
Tel (work) +27 (0)11 717 64 04
Tel (home) +27 (0)18 297 44 51
Fax +27 (0)18 299 24 64
Email [EMAIL PROTECTED]
Web www.oomvanlieshout.net/sander
__
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
Re: [R] Replacing for loop with tapply!?
Thanks Dimitris,
Very impressive! Much faster than before.
Thanks to new found R.basic, I can simply rotate the result with
rotate270{R.basic}:
mat - matrix(sample(-15:50, 365 * 15000, TRUE), 365, 15000)
temps - c(37, 39, 41)
#
#ind - matrix(0, length(temps), ncol(mat))
ind - matrix(0, 4, ncol(mat))
(startDate - date())
[1] Fri Jun 10 12:08:01 2005
for(i in seq(along = temps)) ind[i, ] - colSums(mat temps[i])
ind[4, ] - colMeans(max(mat))
Error in colMeans(max(mat)) : 'x' must be an array of at least two
dimensions
(endDate - date())
[1] Fri Jun 10 12:08:02 2005
ind - rotate270(ind)
ind[1:10,]
V4 V3 V2 V1
1 0 56 75 80
2 0 46 53 60
3 0 50 58 67
4 0 60 72 80
5 0 59 68 76
6 0 55 67 74
7 0 62 77 93
8 0 45 57 67
9 0 57 68 75
10 0 61 66 76
However, I have not managed to get the row maximum using your method? It
should be 50 for most rows, but my first guess code gives an error!
Any suggestions?
Sander
Dimitris Rizopoulos wrote:
maybe you are looking for something along these lines:
mat - matrix(sample(-15:50, 365 * 15000, TRUE), 365, 15000)
temps - c(37, 39, 41)
#
ind - matrix(0, length(temps), ncol(mat))
for(i in seq(along = temps)) ind[i, ] - colSums(mat temps[i])
ind
I hope it helps.
Best,
Dimitris
Dimitris Rizopoulos
Ph.D. Student
Biostatistical Centre
School of Public Health
Catholic University of Leuven
Address: Kapucijnenvoer 35, Leuven, Belgium
Tel: +32/16/336899
Fax: +32/16/337015
Web: http://www.med.kuleuven.ac.be/biostat/
http://www.student.kuleuven.ac.be/~m0390867/dimitris.htm
- Original Message -
From: Sander Oom [EMAIL PROTECTED]
To: r-help@stat.math.ethz.ch
Sent: Friday, June 10, 2005 10:50 AM
Subject: [R] Replacing for loop with tapply!?
Dear all,
We have a large data set with temperature data for weather stations
across the globe (15000 stations).
For each station, we need to calculate the number of days a certain
temperature is exceeded.
So far we used the following S code, where mat88 is a matrix
containing
rows of 365 daily temperatures for each of 15000 weather stations:
m - 37
n - 2
outmat88 - matrix(0, ncol = 4, nrow = nrow(mat88))
for(i in 1:nrow(mat88)) {
# i - 3
row1 - as.data.frame(df88[i, ])
temprow37 - select.rows(row1, row1 m)
temprow39 - select.rows(row1, row1 m + n)
temprow41 - select.rows(row1, row1 m + 2 * n)
outmat88[i, 1] - max(row1, na.rm = T)
outmat88[i, 2] - count.rows(temprow37)
outmat88[i, 3] - count.rows(temprow39)
outmat88[i, 4] - count.rows(temprow41)
}
outmat88
We have transferred the data to a more potent Linux box running R,
but
still hope to speed up the code.
I know a for loop should be avoided when looking for speed. I also
know
the answer is in something like tapply, but my understanding of
these
commands is still to limited to see the solution. Could someone show
me
the way!?
Thanks in advance,
Sander.
--
Dr Sander P. Oom
Animal, Plant and Environmental Sciences,
University of the Witwatersrand
Private Bag 3, Wits 2050, South Africa
Tel (work) +27 (0)11 717 64 04
Tel (home) +27 (0)18 297 44 51
Fax +27 (0)18 299 24 64
Email [EMAIL PROTECTED]
Web www.oomvanlieshout.net/sander
__
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
__
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
--
Dr Sander P. Oom
Animal, Plant and Environmental Sciences,
University of the Witwatersrand
Private Bag 3, Wits 2050, South Africa
Tel (work) +27 (0)11 717 64 04
Tel (home) +27 (0)18 297 44 51
Fax +27 (0)18 299 24 64
Email [EMAIL PROTECTED]
Web www.oomvanlieshout.net/sander
__
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
Re: [R] Replacing for loop with tapply!?
Dear all,
Dimitris and Andy, thanks for your great help. I have progressed to the
following code which runs very fast and effective:
mat - matrix(sample(-15:50, 15 * 10, TRUE), 15, 10)
mat[mat45] - NA
mat-NA
mat
temps - c(35, 37, 39)
ind - rbind(
t(sapply(temps, function(temp)
rowSums(mat temp, na.rm=TRUE) )),
rowSums(!is.na(mat), na.rm=FALSE),
apply(mat, 1, max, na.rm=TRUE))
ind - t(ind)
ind
However, some weather stations have missing values for the whole year.
Unfortunately, the code breaks down (when uncommenting mat-NA).
I have tried 'ifelse' statements in the functions, but it becomes even
more of a mess. I could subset the matrix before hand, but this would
mean merging with a complete matrix afterwards to make it compatible
with other years. That would slow things down.
How can I make the code robust for rows containing all missing values?
Thanks for your help,
Sander.
Dimitris Rizopoulos wrote:
for the maximum you could use something like:
ind[, 1] - apply(mat, 2, max)
I hope it helps.
Best,
Dimitris
Dimitris Rizopoulos
Ph.D. Student
Biostatistical Centre
School of Public Health
Catholic University of Leuven
Address: Kapucijnenvoer 35, Leuven, Belgium
Tel: +32/16/336899
Fax: +32/16/337015
Web: http://www.med.kuleuven.ac.be/biostat/
http://www.student.kuleuven.ac.be/~m0390867/dimitris.htm
- Original Message -
From: Sander Oom [EMAIL PROTECTED]
To: Dimitris Rizopoulos [EMAIL PROTECTED]
Cc: r-help@stat.math.ethz.ch
Sent: Friday, June 10, 2005 12:10 PM
Subject: Re: [R] Replacing for loop with tapply!?
Thanks Dimitris,
Very impressive! Much faster than before.
Thanks to new found R.basic, I can simply rotate the result with
rotate270{R.basic}:
mat - matrix(sample(-15:50, 365 * 15000, TRUE), 365, 15000)
temps - c(37, 39, 41)
#
#ind - matrix(0, length(temps), ncol(mat))
ind - matrix(0, 4, ncol(mat))
(startDate - date())
[1] Fri Jun 10 12:08:01 2005
for(i in seq(along = temps)) ind[i, ] - colSums(mat temps[i])
ind[4, ] - colMeans(max(mat))
Error in colMeans(max(mat)) : 'x' must be an array of at least two
dimensions
(endDate - date())
[1] Fri Jun 10 12:08:02 2005
ind - rotate270(ind)
ind[1:10,]
V4 V3 V2 V1
1 0 56 75 80
2 0 46 53 60
3 0 50 58 67
4 0 60 72 80
5 0 59 68 76
6 0 55 67 74
7 0 62 77 93
8 0 45 57 67
9 0 57 68 75
10 0 61 66 76
However, I have not managed to get the row maximum using your
method? It
should be 50 for most rows, but my first guess code gives an error!
Any suggestions?
Sander
Dimitris Rizopoulos wrote:
maybe you are looking for something along these lines:
mat - matrix(sample(-15:50, 365 * 15000, TRUE), 365, 15000)
temps - c(37, 39, 41)
#
ind - matrix(0, length(temps), ncol(mat))
for(i in seq(along = temps)) ind[i, ] - colSums(mat temps[i])
ind
I hope it helps.
Best,
Dimitris
Dimitris Rizopoulos
Ph.D. Student
Biostatistical Centre
School of Public Health
Catholic University of Leuven
Address: Kapucijnenvoer 35, Leuven, Belgium
Tel: +32/16/336899
Fax: +32/16/337015
Web: http://www.med.kuleuven.ac.be/biostat/
http://www.student.kuleuven.ac.be/~m0390867/dimitris.htm
- Original Message -
From: Sander Oom [EMAIL PROTECTED]
To: r-help@stat.math.ethz.ch
Sent: Friday, June 10, 2005 10:50 AM
Subject: [R] Replacing for loop with tapply!?
Dear all,
We have a large data set with temperature data for weather stations
across the globe (15000 stations).
For each station, we need to calculate the number of days a certain
temperature is exceeded.
So far we used the following S code, where mat88 is a matrix
containing
rows of 365 daily temperatures for each of 15000 weather stations:
m - 37
n - 2
outmat88 - matrix(0, ncol = 4, nrow = nrow(mat88))
for(i in 1:nrow(mat88)) {
# i - 3
row1 - as.data.frame(df88[i, ])
temprow37 - select.rows(row1, row1 m)
temprow39 - select.rows(row1, row1 m + n)
temprow41 - select.rows(row1, row1 m + 2 * n)
outmat88[i, 1] - max(row1, na.rm = T)
outmat88[i, 2] - count.rows(temprow37)
outmat88[i, 3] - count.rows(temprow39)
outmat88[i, 4] - count.rows(temprow41)
}
outmat88
We have transferred the data to a more potent Linux box running R,
but
still hope to speed up the code.
I know a for loop should be avoided when looking for speed. I also
know
the answer is in something like tapply, but my understanding of
these
commands is still to limited to see the solution. Could someone
show
me
the way!?
Thanks in advance,
Sander.
--
__
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
Re: [R] plot3d
Sabine,
It helps us to help you if you tell us what you want to do, preferably
with a code example, and what system/version of R you have: type version.
On my PC I get:
version
_
platform i686-pc-linux-gnu
arch i686
os linux-gnu
system i686, linux-gnu
status
major2
minor1.0
year 2005
month04
day 18
language R
I can not find plot3d on my installation (updated/installed all package
from CRAN yesterday). Quite sure you do not want to start out with a
package that is not on CRAN. Are you sure none of the packages on CRAN
provide the functionality you need?
Cheers,
Sander.
Navarre Sabine wrote:
hello,
to use the function plot3d, i should use the package R.basic!
plot3d {R.basic}
If people know exactly a site to load this package, please give me the URL!
Thanks
Sabine
-
[[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
--
Dr Sander P. Oom
Animal, Plant and Environmental Sciences,
University of the Witwatersrand
Private Bag 3, Wits 2050, South Africa
Tel (work) +27 (0)11 717 64 04
Tel (home) +27 (0)18 297 44 51
Fax +27 (0)18 299 24 64
Email [EMAIL PROTECTED]
Web www.oomvanlieshout.net/sander
__
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
Re: [R] plot3d
Now you made me curious!
Installed the package myself as well!
Maybe you can send an example graph to Romain, to be included in the
graph gallery. Then we can all see what plot3d does!
Thanks,
Sander.
Sander Oom wrote:
Sabine,
It helps us to help you if you tell us what you want to do, preferably
with a code example, and what system/version of R you have: type version.
On my PC I get:
version
_
platform i686-pc-linux-gnu
arch i686
os linux-gnu
system i686, linux-gnu
status
major2
minor1.0
year 2005
month04
day 18
language R
I can not find plot3d on my installation (updated/installed all package
from CRAN yesterday). Quite sure you do not want to start out with a
package that is not on CRAN. Are you sure none of the packages on CRAN
provide the functionality you need?
Cheers,
Sander.
Navarre Sabine wrote:
hello,
to use the function plot3d, i should use the package R.basic!
plot3d {R.basic}
If people know exactly a site to load this package, please give me the
URL!
Thanks
Sabine
-
[[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
--
Dr Sander P. Oom
Animal, Plant and Environmental Sciences,
University of the Witwatersrand
Private Bag 3, Wits 2050, South Africa
Tel (work) +27 (0)11 717 64 04
Tel (home) +27 (0)18 297 44 51
Fax +27 (0)18 299 24 64
Email [EMAIL PROTECTED]
Web www.oomvanlieshout.net/sander
__
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
Re: [R] R Graph Gallery : categorization of the graphs
I agree that a wiki to facilitate submission of graph code could be very effective! Still needs to be well protected against vandalism. Seems a regular backup, to facilitate a clean restore, is the best approach. Romain, would you be willing to set up a wiki within the gallery. Think the wiki and gallery should be close to each other in cyber space! I use DokuWiki privately and it works wonders. Guess the database design will require a bit of thought. For each graph it should be possible to enter multiple categories and sub-categories. then the gallery interface should dynamically build galleries using these categories. Not sure if you would want to create 'meta' categories to facilitate multiple categorization approaches!? We could start with a limited list such as suggested by Chris to populate the gallery in the first instance. Pleased to see this progress!! Sander. Chris Evans wrote: On 6 Jun 2005 at 17:48, Sander Oom wrote: ... much snipped ... The whole point of a gallery is to show something to the user before the user knows what he is looking for. The R help functions currently available are hopeless when you have a picture of a graph in your head without knowing the required commands. ... much snipped ... Belief that good graphics are often as important or more important than inferential tests or even CIs was one of the reasons I've moved over the last 15 years from SPSS (still use it a bit 'cos most colleagues do) through SAS (much better, much better graphics) to S+ (same but more so) to R (same and FLOSS!) The demo graphics and the gallery are wonderful visual arguments for R and also great resources to help us learn. Categories that I think might be useful sometimes might be: describes one variable where that is: dichotomous, categorical (n(categories) 2), polytomous (short), polytomous (many levels), or continuous (might allow something on superimposing different referential distributions describes relationship between two variables where: both are dichotomous or polytomous one is ditto, other is continuous (box violin etc: very useful to see good e.g.s of how to get most appropriate boxplots as it's always possible to get good ones but not always obvious) (pointer to back-to-back histogram in Hmisc here) both are continuous (with and without jitter and weighting blobs) describe relationships between more than two variables... However, the gallery idea is a very powerful one and being able to scroll through and drill down is a useful trick that M$ have, I grudgingly admit, used well so could we mimic their galleries from Excel as someone has suggested and perhaps mimic the drop down graphics picker in S+ (I no longer have access). It's not much help but someone could put up the drop down list for newbies coming from SPSS ... ooh, just opened up my copy (11.0.1) and realised there's a gallery there with the following: bar, line, area, pie, high-low, pareto, control, boxplots, error bar, scatter, histogram, normal P-P, normal Q-Q, sequence, autocorrelations, cross-correlations spectral. Never knew that the was there! The drop down list below that has essentially the same list but with the last three under a sub-heading of time series and ROC curves added. I wonder if someone hosted a wiki for a while at least it would get people contributing code for examples for some of these? The results could transfer to the wonderful graphics gallery as they accumulated. My skills aren't that hot but I'd throw in a few things happily and I'm sure a reward would be hearing of better ways to do things both in terms of coding and in terms of better displays/graphics to use. Cheers all, C -- Dr Sander P. Oom Animal, Plant and Environmental Sciences, University of the Witwatersrand Private Bag 3, Wits 2050, South Africa Tel (work) +27 (0)11 717 64 04 Tel (home) +27 (0)18 297 44 51 Fax +27 (0)18 299 24 64 Email [EMAIL PROTECTED] Web www.oomvanlieshout.net/sander __ 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
Re: [R] A long digression on packages
Oooops, already missed one: 5. search of the R mailing lists: http://maths.newcastle.edu.au/~rking/R/ ad5. never used this before. Think Google also does an excellent job finding these if you start search with R. Sander. Sander Oom wrote: Maybe some of this confusion about search opportunities and pros/cons could be avoided if the search page on CRAN (http://cran.r-project.org/search.html) would be extended to cover all main search tools! Quickly scanning the discussion, I found these: 1- simply Google: some tips and tricks have been mentioned and would be usefully for most users; 2- R site search (external to CRAN) http://finzi.psych.upenn.edu/search.html; 3- from R prompt: help.search(); 4- browser supported search through local help files: R/doc/html/search/SearchEngine.html. ad1. Google is normally my first source using broad keywords for a method or problem. ad2. just discovered this today! ad3. help.search() provides a simple overview, printing the command with the providing package: help.search(rose) Help files with alias or concept or title matching rose using regular expression matching: hirose(boot)Failure Time of PET Film rose.diag(CircStats)Rose Diagram rose.diag(circular) Rose Diagram windrose(circular) Windrose Generator rosavent(climatol) Wind-rose plot Kinship82(clue) Rosenberg-Kim Kinship Terms Partition Data HolidayDatabase(fCalendar) Holiday Calendars and Utilities conc(ineq) Concentration Measures Type 'help(FOO, package = PKG)' to inspect entry 'FOO(PKG) TITLE'. ad4. when installing all package locally, the results produced by the browser supported search can be overwhelming. Even searching within the results often does not help. If commands were printed along with the providing package would be a good improvement. Then the apparent random order of commands listed might also reveal some order. Hope this is a useful addition to the debate, Sander. Chris Evans wrote: On 5 Jun 2005 at 18:44, Jari Oksanen wrote: There are diverse opinions about netiquette. One of the most basic, in my opinion, is this: if someone posts starts a discussion in a certain forum, you shall not divert it to another forum where it may be hidden by most readers, perhaps even by the originator of the thread. With the greatest of respect for Duncan and the R-devel list, I think Jari has a point here. This is one of the most important issues I've seen raised on this list (R-help) in recent months and I think it may be a structural problem for the development of R, in common with that of much FLOSS s'ware, that there's a separation of users and authors that needs thought. There are no perfect answers but too big a separation and projects go techno and it's hard for those of us who can't code C and who are mere users to help those outstanding people on whom we depend hear what we need: sometimes they are so clever, so specialised in their knowledge, or simply in the realm of genius not the ordinary, that they can't see our problem. I have slowly come to respect that a pretty brusque style from our authorities is the only way to prevent this list being a madhouse but I think that Jim's point may fall into that class that is worth some duplicate bandwidth here. I know I've found the problem Jim highlights very confusing and unhelpful at times. Views, which I didn't know, seem helpful but not a real solution to the key problem: they may tangentially help by ensuring that if your needs fit into a view, it becomes more likely that you'll install the packages you need and a local search may tell you what you need. I've taken the inefficient route which suits me of installing just about every package to make it less likely I'll miss something of use to me. That means my search for kappa and Cohen (with ignore.case=FALSE) turns up at least three implementations of aspects of Cohen's kappa. It may already exist, but a web interface that did a help.search over all the packages in the current release version would be great. (If it does exist, sorry, but I'm no dunce and use R nearly every day and try to read much of r-help every day and don't know it, which may say something!) I think there may be a need for some R improvement and automated updating of what I think is Frank Harrell's function finder: http://biostat.mc.vanderbilt.edu/s/finder/finder.html Though I'm not absolutely sure how fitting your works into something like that could be imposed on developers! Another thing that might help would be for a system by which ordinary users would volunteer to pair up with developers for packages and try to suggest adaptations of the help and such like that might make the packages more user friendly. I wouldn't want to do that for the whole of a huge and vital package like MASS or Hmisc (or base or stats!) but I'm up for pairing with a developer on a smaller package
Re: [R] Tabelle zu Matrix
Why ask if you already know the answer. Did you try table()?
Other commands related: reshape {stats} and xtabs {stats}
I learned about this only a couple of days ago on this list!
Meanwhile read the posting guide and send your requests in English.
Good luck,
Sander.
Hansi Weissensteiner wrote:
Hallo!
Ich habe ein Textdokument das folgenden Aufbau aufweisst:
String Integer z.B.
AAA10
BBB15
CCC12
BBB13
AAA11
DDD14
Mein Ziel ist es, die Daten in eine Matrix zu schreiben, ohne dabei mit
Schleifen zu arbeiten. Ist dies möglich? Die entsprechende Matrix sollte dann
wie folgt aufgebaut sein:
AAA BBB CCC DDD
10 15 12 14
11 13 0 0
Ich denke an den Befehl table(), oder gibt es sonstige Befehle?
Danke im Vorraus
MfG
Hansi Weissensteiner
__
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
--
Dr Sander P. Oom
Animal, Plant and Environmental Sciences,
University of the Witwatersrand
Private Bag 3, Wits 2050, South Africa
Tel (work) +27 (0)11 717 64 04
Tel (home) +27 (0)18 297 44 51
Fax +27 (0)18 299 24 64
Email [EMAIL PROTECTED]
Web www.oomvanlieshout.net/sander
__
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
Re: [R] R Graph Gallery : categorization of the graphs
Hi Romain, You have stuck your neck out which is great. You are however not responsible for all the work! Let others join in! Thus a wiki could provide opportunities for other people to contribute! I was suggesting that the wiki be used as a platform for people to submit code examples and for people to discuss these examples! Did not think this would replace the things you need to do to get the data into the gallery!! Maybe Nick would like to kick off the Rgraph wiki pages in the existing Rwiki, following the debates we have had about the categorization!? Cheers, Sander. Romain Francois wrote: Le 07.06.2005 12:36, Sander Oom a écrit : I agree that a wiki to facilitate submission of graph code could be very effective! Still needs to be well protected against vandalism. Seems a regular backup, to facilitate a clean restore, is the best approach. Romain, would you be willing to set up a wiki within the gallery. Think the wiki and gallery should be close to each other in cyber space! Sander, A lot of things are done after getting the code. Syntax highlighting for example. That can't (at least for now) be done online. Certain things will be wiki-like, the WISHLIST can be filled online, in a few days the keywors thing will come and users of the gallery will be able to fill it. Later, people will have the possibility to write a few words on a graph and not just give a note. Nevertheless, I don't think that all should be wikied. I use DokuWiki privately and it works wonders. Guess the database design will require a bit of thought. For each graph it should be possible to enter multiple categories and sub-categories. then the gallery interface should dynamically build galleries using these categories. Not sure if you would want to create 'meta' categories to facilitate multiple categorization approaches!? The keyword thing appears to be really popular and I think it's a good idea. Those keywords may be categorized in the future. That shouldn't be that hard to achieve that. The problem is time. R Graph Gallery is just something I am having fun with. ;-) We could start with a limited list such as suggested by Chris to populate the gallery in the first instance. Pleased to see this progress!! Sander. Thanks to all people who responded to that thread. It really helped me understanding how the things have to be done. Regards, Romain -- Dr Sander P. Oom Animal, Plant and Environmental Sciences, University of the Witwatersrand Private Bag 3, Wits 2050, South Africa Tel (work) +27 (0)11 717 64 04 Tel (home) +27 (0)18 297 44 51 Fax +27 (0)18 299 24 64 Email [EMAIL PROTECTED] Web www.oomvanlieshout.net/sander __ 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
Re: [R] A long digression on packages
Maybe some of this confusion about search opportunities and pros/cons could be avoided if the search page on CRAN (http://cran.r-project.org/search.html) would be extended to cover all main search tools! Quickly scanning the discussion, I found these: 1- simply Google: some tips and tricks have been mentioned and would be usefully for most users; 2- R site search (external to CRAN) http://finzi.psych.upenn.edu/search.html; 3- from R prompt: help.search(); 4- browser supported search through local help files: R/doc/html/search/SearchEngine.html. ad1. Google is normally my first source using broad keywords for a method or problem. ad2. just discovered this today! ad3. help.search() provides a simple overview, printing the command with the providing package: help.search(rose) Help files with alias or concept or title matching rose using regular expression matching: hirose(boot)Failure Time of PET Film rose.diag(CircStats)Rose Diagram rose.diag(circular) Rose Diagram windrose(circular) Windrose Generator rosavent(climatol) Wind-rose plot Kinship82(clue) Rosenberg-Kim Kinship Terms Partition Data HolidayDatabase(fCalendar) Holiday Calendars and Utilities conc(ineq) Concentration Measures Type 'help(FOO, package = PKG)' to inspect entry 'FOO(PKG) TITLE'. ad4. when installing all package locally, the results produced by the browser supported search can be overwhelming. Even searching within the results often does not help. If commands were printed along with the providing package would be a good improvement. Then the apparent random order of commands listed might also reveal some order. Hope this is a useful addition to the debate, Sander. Chris Evans wrote: On 5 Jun 2005 at 18:44, Jari Oksanen wrote: There are diverse opinions about netiquette. One of the most basic, in my opinion, is this: if someone posts starts a discussion in a certain forum, you shall not divert it to another forum where it may be hidden by most readers, perhaps even by the originator of the thread. With the greatest of respect for Duncan and the R-devel list, I think Jari has a point here. This is one of the most important issues I've seen raised on this list (R-help) in recent months and I think it may be a structural problem for the development of R, in common with that of much FLOSS s'ware, that there's a separation of users and authors that needs thought. There are no perfect answers but too big a separation and projects go techno and it's hard for those of us who can't code C and who are mere users to help those outstanding people on whom we depend hear what we need: sometimes they are so clever, so specialised in their knowledge, or simply in the realm of genius not the ordinary, that they can't see our problem. I have slowly come to respect that a pretty brusque style from our authorities is the only way to prevent this list being a madhouse but I think that Jim's point may fall into that class that is worth some duplicate bandwidth here. I know I've found the problem Jim highlights very confusing and unhelpful at times. Views, which I didn't know, seem helpful but not a real solution to the key problem: they may tangentially help by ensuring that if your needs fit into a view, it becomes more likely that you'll install the packages you need and a local search may tell you what you need. I've taken the inefficient route which suits me of installing just about every package to make it less likely I'll miss something of use to me. That means my search for kappa and Cohen (with ignore.case=FALSE) turns up at least three implementations of aspects of Cohen's kappa. It may already exist, but a web interface that did a help.search over all the packages in the current release version would be great. (If it does exist, sorry, but I'm no dunce and use R nearly every day and try to read much of r-help every day and don't know it, which may say something!) I think there may be a need for some R improvement and automated updating of what I think is Frank Harrell's function finder: http://biostat.mc.vanderbilt.edu/s/finder/finder.html Though I'm not absolutely sure how fitting your works into something like that could be imposed on developers! Another thing that might help would be for a system by which ordinary users would volunteer to pair up with developers for packages and try to suggest adaptations of the help and such like that might make the packages more user friendly. I wouldn't want to do that for the whole of a huge and vital package like MASS or Hmisc (or base or stats!) but I'm up for pairing with a developer on a smaller package if anyone thinks that would be helpful. Thoughts for what they're worth. Thanks a million to all developers ... asbestos suit on! Chris -- Dr Sander P. Oom Animal, Plant and Environmental
Re: [R] Polar Graph
Navarre Sabine wrote:
Hi,
I would like to do a polar graph (=star graph) ! is that graph existing on R?
Because more softwares can do that but I don't found it on R!
Thanks
Sabine
-
ils, photos et vids !
[[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
Try: rose.diag {CircStats}!
Sander.
--
Dr Sander P. Oom
Animal, Plant and Environmental Sciences,
University of the Witwatersrand
Private Bag 3, Wits 2050, South Africa
Tel (work) +27 (0)11 717 64 04
Tel (home) +27 (0)18 297 44 51
Fax +27 (0)18 299 24 64
Email [EMAIL PROTECTED]
Web www.oomvanlieshout.net/sander
__
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
Re: [R] R Graph Gallery : categorization of the graphs
Barry Rowlingson wrote: Romain Francois wrote: Graphics will be classified in : - categories - sub-categories within those categories So far so good. Maybe, maybe not! Would a system of keywords work better than strict hierarchical categories, as long as plots can have more than one keyword attached? Someone might be interested in all 3d plots, or all plots related to model fitting, and these categories may not be distinct! Barry __ 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 The whole point of a gallery is to show something to the user before the user knows what he is looking for. The R help functions currently available are hopeless when you have a picture of a graph in your head without knowing the required commands. Thus a gallery based on some categorization is required. The gallery already has a search facility for 'un-categorized' searching. At the same time the gallery will provide cross links between different graphs using the same command/package or other similarities. I suggested that the categorization for sigmaplot graphs could be a useful start: http://www.systat.com/products/SigmaPlot/productinfo/?sec=1003 The gallery is based on a database, so keyword searching is no problem. Of course graphs will likely be members of more then one sub-category. This is also no problem thanks to the database approach! Hopefully we can agree on a categorization, such that we can start populating the gallery! Cheers, Sander. -- Dr Sander P. Oom Animal, Plant and Environmental Sciences, University of the Witwatersrand Private Bag 3, Wits 2050, South Africa Tel (work) +27 (0)11 717 64 04 Tel (home) +27 (0)18 297 44 51 Fax +27 (0)18 299 24 64 Email [EMAIL PROTECTED] Web www.oomvanlieshout.net/sander __ 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
Re: [R] How to 'de-cross' a table?
Hi Thomas, Your code works perfectly! Thanks a lot, Sander. Thomas Lumley wrote: On Fri, 3 Jun 2005, Sander Oom wrote: Dear R users, I have received a table in the following format: id a b c1 c2 d1 d2 1 1 1 65 97 78 98 2 1 2 65 97 42 97 3 2 1 65 68 97 98 4 2 2 65 97 97 98 Factors of the design are: a, b, and e, where e has levels c and d. The levels c and d then have 2 replicates (r) each. Now I would like to get: id a b e r value 1 1 1 c 1 65 2 1 1 c 2 97 3 1 1 d 1 78 4 1 1 d 2 98 Any suggestions on how to tackle this? I'm not sure what the 'task' is called, so struggle to effectively search the web or R help. reshape() is the probably function. It seems you need to use it twice, for the two levels of uncrossing d1-reshape(d, direction=long, varying=list(c(c1,d1),c(c2,d2)), v.names=c(rep1,rep2),timevar=r, times=c(c,d)) d1 id a b r rep1 rep2 1.c 1 1 1 c 65 97 2.c 2 1 2 c 65 97 3.c 3 2 1 c 65 68 4.c 4 2 2 c 65 97 1.d 1 1 1 d 78 98 2.d 2 1 2 d 42 97 3.d 3 2 1 d 97 98 4.d 4 2 2 d 97 98 reshape(d1, direction=long, varying=list(c(rep1,rep2)), v.names=value,idvar=tmp, timevar=r) id a b r value tmp 1.1 1 1 1 165 1 2.1 2 1 2 165 2 3.1 3 2 1 165 3 4.1 4 2 2 165 4 5.1 1 1 1 178 5 6.1 2 1 2 142 6 7.1 3 2 1 197 7 8.1 4 2 2 197 8 1.2 1 1 1 297 1 2.2 2 1 2 297 2 3.2 3 2 1 268 3 4.2 4 2 2 297 4 5.2 1 1 1 298 5 6.2 2 1 2 297 6 7.2 3 2 1 298 7 8.2 4 2 2 298 8 You can then delete the `tmp` variable if you don't need it. An easier way to uncross would be with stack() stack(d[,4:7]) values ind 1 65 c1 2 65 c1 3 65 c1 4 65 c1 5 97 c2 6 97 c2 7 68 c2 8 97 c2 9 78 d1 10 42 d1 11 97 d1 12 97 d1 13 98 d2 14 97 d2 15 98 d2 16 98 d2 but you lose the other variables. -thomas __ 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 -- Dr Sander P. Oom Animal, Plant and Environmental Sciences, University of the Witwatersrand Private Bag 3, Wits 2050, South Africa Tel (work) +27 (0)11 717 64 04 Tel (home) +27 (0)18 297 44 51 Fax +27 (0)18 299 24 64 Email [EMAIL PROTECTED] Web www.oomvanlieshout.net/sander __ 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
[R] How to 'de-cross' a table?
Dear R users, I have received a table in the following format: id a b c1 c2 d1 d2 1 1 1 65 97 78 98 2 1 2 65 97 42 97 3 2 1 65 68 97 98 4 2 2 65 97 97 98 Factors of the design are: a, b, and e, where e has levels c and d. The levels c and d then have 2 replicates (r) each. Now I would like to get: id a b e r value 1 1 1 c 1 65 2 1 1 c 2 97 3 1 1 d 1 78 4 1 1 d 2 98 Any suggestions on how to tackle this? I'm not sure what the 'task' is called, so struggle to effectively search the web or R help. Thanks, Sander. Dr Sander P. Oom Animal, Plant and Environmental Sciences, University of the Witwatersrand Private Bag 3, Wits 2050, South Africa Tel (work) +27 (0)11 717 64 04 Tel (home) +27 (0)18 297 44 51 Fax +27 (0)18 299 24 64 Email [EMAIL PROTECTED] Web www.oomvanlieshout.net/sander __ 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
Re: [R] Increasing Console Paste Buffer
An interesting thought just came to me when reading this discussion! I use both R and Latex and have never had the trouble of overlooking error messages when debugging long Latex code! Of course this is because when compiling a latex document, a summary of the compilation process is provided at the end! If any errors occurred, they will be mentioned in the summary. Maybe R could provide the same summary as an optional part of the source() command!? Cheers, Sander. Gavin Simpson wrote: Jan T. Kim wrote: On Tue, May 31, 2005 at 11:47:05PM +0100, Gavin Simpson wrote: Manuel Morales wrote: Hello list. I'm using R from the gnome-terminal in Fedora. My preference is to write programs in VIM, and then source the file from R, or copy and paste the lines into the console. I'm wondering if there is a way to increase the paste buffer as an alternative to sourcing large analyses. As was mentioned in a recent thread on Linux GUI's, I find that if I paste in a large amount of text, the lines end up getting cut off at some point. I wonder if this is an R restriction, because it seems like I am able to paste substantially more text in other console-based programs. Is there any way to increase the amount of text that I can paste into an R session? Thanks! Manuel Manuel, Maybe I misunderstand what you mean by lines end up getting cut off at some point so correct me if I got it wrong, but I assume you mean that after a certain number of lines entered you can no longer scroll back up and view the earlier lines? I think that this is not an issue of the scroll buffer, but of buffers internal to the terminal program or the shell, which are designed to hold keyboard input and which can be overwhelmed by the rate of input when large text selections are pasted in, as this appears as though thousands of keys had been typed almost instantaneously from their view, so to speak. I did say I was guessing :-) For these reasons, I generally strongly recommend against pasting into terminals. Thanks for this Jan. I haven't noticed this myself but then again I hate copy/paste and rarely use R outside emacs/ess these days. In R, use the source() instead... ;-) Agreed. source(filename, echo = TRUE) will sort of replicate the behaviour the original poster would get if they like to see the commands printed among the results. But if he is pasting in that much data, Manuel will still have to increase the buffer on the terminal, especially if he is using one of the defaults in FC3 as the output will quickly get lost. Best regards, Jan All the best, Gav -- Dr Sander P. Oom Animal, Plant and Environmental Sciences, University of the Witwatersrand Private Bag 3, Wits 2050, South Africa Tel (work) +27 (0)11 717 64 04 Tel (home) +27 (0)18 297 44 51 Fax +27 (0)18 299 24 64 Email [EMAIL PROTECTED] Web www.oomvanlieshout.net/sander __ 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
Re: [R] Increasing Console Paste Buffer
Indeed it does! Sorry for the impulsive response!
Sander.
Duncan Murdoch wrote:
Sander Oom wrote:
An interesting thought just came to me when reading this discussion! I
use both R and Latex and have never had the trouble of overlooking
error messages when debugging long Latex code!
Of course this is because when compiling a latex document, a summary
of the compilation process is provided at the end! If any errors
occurred, they will be mentioned in the summary.
Maybe R could provide the same summary as an optional part of the
source() command!?
I think it does, doesn't it? R will stop at the first error and print
it, e.g.
source('c:/temp/test.R')
Error in parse(file, n = -1, NULL, ?) : syntax error on line 4
If there were only warnings, it will show them at the end:
source('c:/temp/test.R')
Warning messages:
1: longer object length
is not a multiple of shorter object length in: 1:3 + 1:4
2: longer object length
is not a multiple of shorter object length in: 1:3 + 1:4
Even if you use echo=TRUE, these summaries show up at the end. It's
only if you use cut and paste that you might miss these.
Duncan Murdoch
__
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
--
Dr Sander P. Oom
Animal, Plant and Environmental Sciences,
University of the Witwatersrand
Private Bag 3, Wits 2050, South Africa
Tel (work) +27 (0)11 717 64 04
Tel (home) +27 (0)18 297 44 51
Fax +27 (0)18 299 24 64
Email [EMAIL PROTECTED]
Web www.oomvanlieshout.net/sander
__
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
Re: [R] R GUI for Linux?
Hi Charles, Warm felt sympathies for your struggles. I consider myself a happy GUI user and have also struggled with the 'command line' history and lack of out-of-the-box functionality associated with Linux. However, Linux does have many, many advantages over other OS's, so I will stick to it. That said, I am reluctantly suggesting it to others, as there are lots of other things one can do on a Saturday afternoon. ;-) After many long afternoons, I have come to the conclusion that accepting what comes as standard is the best approach to using Linux if you do have better things to do. I tried ESS, but I found it impenetrable at my first try, so gave up. I run SuSE linux and R without any problem. The RPM was downloaded from CRAN and installed without any errors. I run update.packages to download, install and refresh all contributed and other packages. I also looked for an R GUI, but must admit I stopped way before you did. I ended up using my favorite editor (Kate, comes standard with KDE) and cut and paste code into an R console (available as standard in the Kate window). Gives a nice clean window for writing code and running R. Of course Kate has code highlighting facilities (as standard). I soon realized that the console has a limited buffer for commands, such that long code sequences are abruptly ended part way through. Thus I have a second console open from which I source whole R files. The setup gives a relatively comfortable code debugging environment! I was amazed when I read the instruction for JGR on Linux. I thought the whole point of Java was to create platform independent software. I have given up on any instructions that tell me to run 'make'. Linux distributions are just to idiosyncratic for it to be worth the effort. I'll just wait until Novell packs JGR with the new version of SuSE Linux. Even if I have to pay something for the added bonus! Will keep following the GUI developments with interest, Sander. White, Charles E WRAIR-Wash DC wrote: I feel your pain. grin I am a new Linux user who has spent most of the weekend trying to get a functional R setup. When I installed Fedora Core 3 (FC3) on my home computer, I thought using R in a terminal would be a snap. I installed R using the rpm packages and tried to use it with the FC3 default terminal (GNOME Terminal 2.7.3). Before long, I found out the terminal was rudely discarding output beyond a set number of lines. I could increase the number of lines kept by the terminal but that didn't strike me as an acceptable solution. Cutting to my stress relieving intermediate solution, I am currently using xemacs with ESS as my R programming environment under FC3. Eventually, I will want to run JGR as my programming environment and Rcmdr as both a teaching tool and means to distribute code to some of my clients. On my way to xemacs, I also tried to install emacs and gnomeGUI. I will briefly document my experience with trying to install each of these packages below: XEMACS with ESS: XEMACS is within the 'walled garden' of packages tuned specifically to run under FC3 and XEMACS provides a tuned installation for ESS. Since I had already compiled R from source with shared libraries enabled (the rmp does not enable shared libraries), I don't know if XEMACS will work with the rpm version of R. Note also that I installed this package using yum; 'Add or Remove Applications' lists xemacs but wouldn't allow me to install. JGR: I have installed jdk1.5.0_03 and MOST of the the output from make looks like JGR is compiling correctly. JavaGD and rJava are not finding jini.h. I don't see an explicit statement of how to start JGR but I assume that is done by typing JGR (or maybe jgr) in a terminal window. Nothing happens. Two potential problems are: (a) I never should have downloaded JavaGD and rJava from CRAN (they won't uninstall, deleting the directories doesn't stop the problem, and I can't use yum to remove R to start over because yum doesn't recognise that R is installed.) or (b) I need to uninstall some of the stray versions of java littering my hard drive. I haven't removed the symbolic link between jre1.5.0_02 and firefox. Rcmdr: There are all sorts of things in FC3 that seem to be tcl/tk related but Rcmdr doesn't seem to work with them. Since some are part of the base FC3 installation, I'm nervious about replacing them or installing competing software. Potentially conflicting software in FC3 are listed below: tcl.i386 8.4.7-2installed tclx.i3868.3.5-4installed db4-tcl.i386 4.2.52-6 base postgresql-tcl.i386 7.4.8-1.FC3.1 updates-released ruby-tcltk.i386 1.8.2-1.FC3.1 updates-released tcl-devel.i386 8.4.7-2base tcl-html.i386
Re: [R] R GUI for Linux?
HI Robert, Of course Linux already has a console! Just type R in the Terminal console and R will start (assuming all is installed correctly). Graphics will be launched in separate windows. If you want more then the Terminal console, try: http://www.sciviews.org/_rgui/ Good luck, Sander. Robert Citek wrote: Hello all, I noticed that both Windows and OS X version of R have a GUI (Rconsole). Is there a GUI for Linux? I'm running Debian on which the CLI for R works just fine. Regards, - Robert __ 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 -- Dr Sander P. Oom Animal, Plant and Environmental Sciences, University of the Witwatersrand Private Bag 3, Wits 2050, South Africa Tel (work) +27 (0)11 717 64 04 Tel (home) +27 (0)18 297 44 51 Fax +27 (0)18 299 24 64 Email [EMAIL PROTECTED] Web www.oomvanlieshout.net/sander __ 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
Re: [R] Soil texture triangle in R?
Hi Jim,
Your email was classified as spam, so I missed it previously.
Unfortunately you did not send a cc to the R-help list. I send a cc to
the mailing list now, so all code gets archived.
Thanks for the improvements on your function! It seems that drawing the
ternary graph and the points with the generic plot functions works well.
Makes the function less dependent on other packages!
Will merge the two functions into one and post it back to the mailing
list! Then the graph might be ready for the graph gallery!
Thanks,
Sander.
Jim Lemon wrote:
Sander Oom wrote:
Hi Jim,
This looks impressive! It gives me the 'background' graph. However,
I'm not sure how I can use this function to plot my soil texture
values! Can you explain?
I would like to be able to plot my soil texture samples in the same
graph as the one your function plots.
Yes, that's just the background figure for which you attached the link.
I was unsure of how you were going to use this, that is, do you just
place a symbol for each soil sample? Aha! Just read your latest email
and I think that's what you want.
Let's say you have one or more soil samples with the proportions of
clay, sand and silt.
sample1-c(0.1,0.4,0.5)
sample2-c(0.2,0.5,0.3)
soil.samples-rbind(sample1,sample2)
colnames(soil.samples)-c(clay,sand,silt)
This more complicated function allows you to overlay colored points
representing soil samples on either an empty triangle, a gridded
triangle or a triangle with soil type names. It also has an optional
legend.
par(ask=TRUE)
soil.triangle(soil.samples)
soil.triangle(soil.samples,show.grid=TRUE)
soil.triangle(soil.samples,soil.names=TRUE,legend=TRUE)
par(ask=FALSE)
Jim
soil.triangle-function(soilprop,pch=NULL,col=NULL,soil.names=FALSE,
show.grid=FALSE,show.legend=FALSE) {
if(missing(soilprop))
stop(Usage:
soil.triangle(soilprop,pch=NULL,col=NULL,soil.names=FALSE,show.grid=FALSE))
if(!is.matrix(soilprop))
stop(soilprop must be a matrix with at least three columns and one
row.)
if(any(soilprop 1) || any(soilprop 0))
stop(All soil proportions must be between zero and one.)
if(any(abs(rowSums(soilprop)-1) 0.01))
warning(At least one set of soil proportions does not equal one.)
oldpar-par(no.readonly=TRUE)
plot(0:1,type=n,axes=FALSE,xlim=c(0,1.1),ylim=c(0,1),
main=Soil Triangle,xlab=,ylab=)
# first draw the triangle
x1-c(0,0,0.5)
sin60-sin(pi/3)
x2-c(1,0.5,1)
y1-c(0,0,sin60)
y2-c(0,sin60,0)
segments(x1,y1,x2,y2)
# now the bottom internal ticks
bx1-seq(0.1,0.9,by=0.1)
bx2-bx1-0.01
by1-rep(0,9)
by2-rep(0.02*sin60,9)
segments(bx1,by1,bx2,by2)
text(bx1,by1-0.03,as.character(rev(seq(10,90,by=10
# now the left internal ticks
ly1-bx1*sin60
lx1-bx1*0.5
lx2-lx1+0.02
ly2-ly1
segments(lx1,ly1,lx2,ly2)
text(lx1-0.03,ly1,as.character(seq(10,90,by=10)))
# right internal ticks
rx1-rev(lx1+0.5-0.01)
rx2-rx1+0.01
ry1-ly1-0.02*sin60
ry2-ly2
segments(rx1,ry1,rx2,ry2)
if(show.grid) {
segments(bx2,by2,lx1,ly1,lty=3)
segments(lx2,ly2,rx2,ry2,lty=3)
segments(rev(rx1),rev(ry1),bx1,by1,lty=3)
}
text(rx2+0.03,ry1+0.025,as.character(rev(seq(10,90,by=10
text(0.5,0.9,100% clay)
par(xpd=TRUE)
text(-0.1,0,100% sand)
text(1.1,0,100% loam)
text(0.07,0.43,percent clay)
text(0.93,0.43,percent silt)
text(0.5,-0.1,percent sand)
if(soil.names) {
# boundary of clay with extensions
x1-c(0.275,0.35,0.6)
x2-c(0.4,0.79,0.7)
y1-c(0.55*sin60,0.41*sin60,0.41*sin60)
y2-c(0.285*sin60,0.41*sin60,0.6*sin60)
segments(x1,y1,x2,y2)
# lower bound of clay loam silty divider
x1-c(0.4,0.68)
x2-c(0.86,0.6)
y1-c(0.285*sin60,0.285*sin60)
y2-c(0.285*sin60,0.41*sin60)
segments(x1,y1,x2,y2)
x1-c(0.185,0.1,0.37)
x2-c(0.36,0.37,0.4)
y1-c(0.37*sin60,0.2*sin60,0.2*sin60)
y2-c(0.37*sin60,0.2*sin60,0.285*sin60)
segments(x1,y1,x2,y2)
# sand corner
x1-c(0.05,0.075)
x2-c(0.12,0.3)
y1-c(0.1*sin60,0.15*sin60)
y2-c(0,0)
segments(x1,y1,x2,y2)
x1-c(0.37,0.42,0.5,0.8,0.86)
x2-c(0.42,0.54,0.65,0.86,0.94)
y1-c(0.2*sin60,0.08*sin60,0,0,0.12*sin60)
y2-c(0.08*sin60,0.08*sin60,0.285*sin60,0.12*sin60,0.12*sin60)
segments(x1,y1,x2,y2)
text(0.5,0.57,Clay)
text(0.7,0.49*sin60,Silty)
text(0.7,0.44*sin60,clay)
text(0.73,0.37*sin60,Silty clay)
text(0.73,0.33*sin60,loam)
text(0.5,0.35*sin60,Clay loam)
text(0.27,0.43*sin60,Sandy)
text(0.27,0.39*sin60,clay)
text(0.27,0.3*sin60,Sandy clay)
text(0.27,0.26*sin60,loam)
text(0.25,0.13*sin60,Sandy loam)
text(0.13,0.075*sin60,Loamy)
text(0.15,0.035*sin60,sand)
text(0.055,0.021,Sand)
text(0.49,0.18*sin60,Loam)
text(0.72,0.15*sin60,Silt loam)
text(0.9,0.06*sin60,Silt)
}
if(is.null(pch)) pch-1:nrow(soilprop)
if(is.null(col)) col-2:(nrow(soilprop)+1)
points(1-soilprop[,2]+(soilprop[,2]-(1-soilprop[,3
Re: [R] Soil texture triangle in R?
Dear R users,
Please find attached a new plot function, plot.soiltexture, to plot soil
texture data on a triangular plot with an optional backdrop of the USDA
soil texture classification, written by Jim Lemon and me.
I tried to write the function and documentation confirm the R
conventions. However, this is a new experience for me, so any comments
and suggestions are welcome!
I tried to find a suitable package for the plot function, but none are
obvious.
I have approached Todd Skaggs to ask permission to include sample data
from the paper:
Skaggs, T.H., L.M. Arya, P.J. Shouse, and B.P. Mohanty. 2001. Estimating
particle-size distribution from limited soil texture data. Soil Sci.
Soc. Am. J., 65:1038-1044.
http://soil.scijournals.org/cgi/content/full/65/4/1038
Things still to do:
- rotate axis labels;
- rotate axis tick labels
- provide option to plot ticks inside or outside plot area.
Thus making it look like:
http://soils.usda.gov/technical/manual/images/fig3-16_large.jpg
or
http://soil.scijournals.org/content/vol65/issue4/images/large/1038f2.jpeg
Enjoy,
Sander.
# Description
#
# Plots soil texture data on an equilateral triangle. Provides the option to
# draw the USDA soil texture classification as a backdrop.
#
# Usage
#
# plot.soiltexture - function(soiltexture, main=NULL, pch=1, col=black,
# soil.names=TRUE, soil.lines=TRUE, show.points=TRUE,
# show.clabels=FALSE, show.grid=FALSE, show.legend=FALSE,
# col.names=grey, col.lines=grey, bg.pch=white,
# col.grid=grey, lty.grid=3)
#
# Arguments
#
# soiltexture a matrix with three columns (see details)
# main title of the plot; defaults to NULL
# pch variable or vector containing point symbols; defaults to 1
# col variable or vector containing point colour; defaults to black
# soil.namesa logical value indicating whether the soil texture class names
are printed; defaults to TRUE
# soil.linesa logical value indicating whether the soil texture class
division lines are plotted; defaults to TRUE
# show.points a logical value indicating whether the soil sample points are
plotted; defaults to TRUE
# show.clabels a logical value indicating whether the corner labels for 100%
sand, silt and clay are printed; defaults to FALSE
# show.grid a logical value indicating whether the fixed grid is plotted;
defaults to FALSE
# show.legend a logical value indicating whether the legendis plotted;
defaults to FALSE
# col.names colour of the soil texture class names; defaults to grey
# col.lines colour of the soil texture class division lines; defaults to
grey
# bg.pchcolour of the points symbols when pch 21:25 are used; defaults
to white
# col.grid colour of the fixed grid; defaults to grey
# lty.grid line style for the fixed grid; defaults to 3
#
# Details
#
# The object soiltexture must be a matrix with at least three columns and one
row.
# Columns should contain data for Sand, Silt, and Clay, in that order! Sand,
Silt,
# and Clay should be expressed in proportions between 0 and 1.
#
# The soil sample points' coordinates are calculated using simple trigonometry.
# Thus, the coordinates of a point P(Sand,Silt,Clay), where Sand + Silt + Clay
= 1,
# are: P(1-Sand+(Sand-(1-Silt))*0.5, Clay*sin(pi/3))
#
# Author(s)
#
# Jim Lemon
# Sander Oom
#
# References
#
# Soil Survey Division Staff. 1993. Soil survey manual. Soil Conservation
Service.
# U.S. Department of Agriculture Handbook 18.
#
# Examples
#
# # some triangular data
# library(MASS)
# tmp-(Skye/100)
# # colnames choosen to be consistent with MASS-fig4.4
# colnames(tmp) - c(Clay,Sand,Silt)
# soiltexture - cbind(tmp$Sand,tmp$Silt,tmp$Clay)
# # the USDA backdrop in black
# plot.soiltexture(NULL, show.points=FALSE, col.names=black,
col.lines=black)
# # the USDA backdrop and a fixed grid in grey
# plot.soiltexture(NULL, show.points=FALSE, show.grid=TRUE)
# # soil sample points with backdrop in grey
# plot.soiltexture(soiltexture)
plot.soiltexture - function(soiltexture, main=Soil Texture Plot, pch=1,
col=black,
soil.names=TRUE, soil.lines=TRUE,
show.points=TRUE, show.clabels=FALSE, show.grid=FALSE, show.legend=FALSE,
col.names=grey, col.lines=grey, bg.pch=white,
col.grid=grey, lty.grid=3) {
if(show.points) {
## error checking
if(missing(soiltexture))
stop(Usage: plot.soiltexture(soiltexture, pch=NULL, col=NULL,
soil.names=FALSE, show.grid=FALSE))
if(!is.matrix(soiltexture))
stop(Object soiltexture must be a matrix with at least three columns
(for Sand, Silt, and Clay) and one row.)
if(any(soiltexture 1) || any(soiltexture 0))
stop(All soil texture proportions must be between zero and one.)
if(any(abs(rowSums(soiltexture)-1) 0.01))
warning(At least one set of soil texture proportions does not equal
one.)
}
oldpar-par(no.readonly=TRUE)
sin60-sin(pi/3)
par(xpd=TRUE)
# create empty canvas to draw
plot(0:1,type=n,axes=FALSE
[R] Soil texture triangle in R?
Dear R users, has anybody made an attempt to create the soil texture triangle graph in R? For an example see here: http://www.teachingkate.org/images/soiltria.gif I would like to get the lines in black and texture labels in gray to allow for plotting my texture results on top. Any examples or suggestions are very welcome! Thanks in advance, Sander. -- Dr Sander P. Oom Animal, Plant and Environmental Sciences, University of the Witwatersrand Private Bag 3, Wits 2050, South Africa Tel (work) +27 (0)11 717 64 04 Tel (home) +27 (0)18 297 44 51 Fax +27 (0)18 299 24 64 Email [EMAIL PROTECTED] Web www.oomvanlieshout.net/sander __ 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
Re: [R] Soil texture triangle in R?
Right,
Got the data points plotted on top of the soil texture background,
thanks to Jim and ternaryplot{vcd}! See code below.
Now there is some fine tuning to do, as it should really look like this
graph:
http://soil.scijournals.org/content/vol65/issue4/images/large/1038f2.jpeg
Things to do:
- rotate axis labels;
- correct small errors in class divisions;
- correct the partial covering of the bottom tick labels;
- rotate ticks in order to simplify viewing the graph.
Any help still appreciated!
Cheers,
Sander.
soil.triangle - function() {
oldpar - par(no.readonly=TRUE)
## now the bottom internal ticks
x1-seq(0.1,0.9,by=0.1)
x2-x1
y1-rep(0,9)
y2-rep(-0.02,9)
segments(x1,y1,x2,y2)
text(x1,y1-0.03,as.character(rev(seq(10,90,by=10 #, cex=0.8)
## now the left internal ticks
y1-x1*sin60
x1-x1*0.5
x2-x1+0.02*sin60
y2-y1-0.02*0.5
segments(x1,y1,x2,y2)
text(x1-0.03,y1+0.015,as.character(seq(10,90,by=10)))
## now the right internal ticks
x1-rev(x1+0.5-0.02*sin60)
x2-x1+0.02*sin60
segments(x1,y2,x2,y1)
text(x2+0.03,y1+0.015,as.character(rev(seq(10,90,by=10
## the labels at the corners
par(xpd=TRUE)
# text(0.5,0.9,100% clay)
# text(-0.1,0,100% sand)
# text(1.1,0,100% loam)
text(0.09,0.43,% Clay)
text(0.90,0.43,% Silt)
text(0.5,-0.1,% Sand)
# boundary of clay with extensions
x1-c(0.275,0.35,0.6)
x2-c(0.4,0.79,0.7)
y1-c(0.55*sin60,0.41*sin60,0.41*sin60)
y2-c(0.285*sin60,0.41*sin60,0.6*sin60)
segments(x1,y1,x2,y2, col=grey)
text(0.5,0.57,Clay, col=grey)
# lower bound of clay loam silty divider
x1-c(0.4,0.68)
x2-c(0.86,0.6)
y1-c(0.285*sin60,0.285*sin60)
y2-c(0.285*sin60,0.41*sin60)
segments(x1,y1,x2,y2, col=grey)
text(0.7,0.49*sin60,Silty, col=grey)
text(0.7,0.44*sin60,clay, col=grey)
text(0.73,0.37*sin60,Silty clay, col=grey)
text(0.73,0.33*sin60,loam, col=grey)
text(0.5,0.35*sin60,Clay loam, col=grey)
x1-c(0.185,0.1,0.37)
x2-c(0.36,0.37,0.4)
y1-c(0.37*sin60,0.2*sin60,0.2*sin60)
y2-c(0.37*sin60,0.2*sin60,0.285*sin60)
segments(x1,y1,x2,y2, col=grey)
text(0.27,0.43*sin60,Sandy, col=grey)
text(0.27,0.39*sin60,clay, col=grey)
text(0.27,0.3*sin60,Sandy clay, col=grey)
text(0.27,0.26*sin60,loam, col=grey)
# sand corner
x1-c(0.05,0.075)
x2-c(0.12,0.3)
y1-c(0.1*sin60,0.15*sin60)
y2-c(0,0)
segments(x1,y1,x2,y2, col=grey)
text(0.25,0.13*sin60,Sandy loam, col=grey)
text(0.13,0.075*sin60,Loamy, col=grey)
text(0.15,0.035*sin60,sand, col=grey)
text(0.055,0.021,Sand, col=grey)
x1-c(0.37,0.42,0.5,0.8,0.86)
x2-c(0.42,0.54,0.65,0.86,0.94)
y1-c(0.2*sin60,0.08*sin60,0,0,0.12*sin60)
y2-c(0.08*sin60,0.08*sin60,0.285*sin60,0.12*sin60,0.12*sin60)
segments(x1,y1,x2,y2, col=grey)
text(0.49,0.18*sin60,Loam, col=grey)
text(0.72,0.15*sin60,Silt loam, col=grey)
text(0.9,0.06*sin60,Silt, col=grey)
par(oldpar)
}
tmp - array(dim=c(10,3))
tmp[,1] - abs(rnorm(10)*20)
tmp[,2] - abs(rnorm(10)*10)
tmp[,3] - 100-tmp[,1]-tmp[,2]
tmp
library(vcd)
## Mark groups
ternaryplot(tmp,
grid=FALSE,
dimnames.position = none,
pch=1, col=black,
scale=1, main=NULL,
prop.size=FALSE,
)
soil.triangle()
Sander Oom wrote:
Hi Jim,
This looks impressive! It gives me the 'background' graph. However, I'm
not sure how I can use this function to plot my soil texture values! Can
you explain?
I would like to be able to plot my soil texture samples in the same
graph as the one your function plots.
Maybe I should try to figure out how to replicate your code inside a
ternaryplot{vcd} call.
Cheers,
Sander.
Jim Lemon wrote:
Sander Oom wrote:
Dear R users,
has anybody made an attempt to create the soil texture triangle graph
in R? For an example see here:
http://www.teachingkate.org/images/soiltria.gif
I would like to get the lines in black and texture labels in gray to
allow for plotting my texture results on top.
Any examples or suggestions are very welcome!
It's not too hard to write a plot function to do this, but I'm not sure
whether this will be what you want. Anyway, try it out.
Jim
soil.triangle-function() {
oldpar-par(no.readonly=TRUE)
plot(0:1,type=n,axes=FALSE,xlim=c(0,1.1),ylim=c(0,1),
main=Soil Triangle,xlab=,ylab=)
# first draw the triangle
x1-c(0,0,0.5)
sin60-sin(pi/3)
x2-c(1,0.5,1)
y1-c(0,0,sin60)
y2-c(0,sin60,0)
segments(x1,y1,x2,y2)
# now the bottom internal ticks
x1-seq(0.1,0.9,by=0.1)
x2-x1
y1-rep(0,9)
y2-rep(0.02,9)
segments(x1,y1,x2,y2)
text(x1,y1-0.03,as.character(rev(seq(10,90,by=10
# now the left internal ticks
y1-x1*sin60
x1-x1*0.5
x2-x1+0.02*sin60
y2-y1-0.02*0.5
segments(x1,y1,x2,y2)
text(x1-0.03,y1+0.015,as.character(seq(10,90,by=10)))
x1-rev(x1+0.5-0.02*sin60)
x2-x1+0.02*sin60
segments(x1,y2,x2,y1)
text(x2+0.03,y1+0.015,as.character(rev(seq(10,90,by=10
text
Re: [R] Soil texture triangle in R?
Cleaned up the class divisions and created a full function.
Still to do:
- rotate axis labels;
- correct the partial covering of the bottom tick labels;
- rotate ticks in order to simplify viewing the graph.
See:
http://soil.scijournals.org/content/vol65/issue4/images/large/1038f2.jpeg
Wonder whether triangle.plot{ade4} will give more flexibility!?
Anyway, hopefully the result so far is useful for other people.
Cheers,
Sander.
plot.soiltexture - function(x,pch,col) {
## triangle plots:
## triangle.plot {ade4}
## triplot{klaR}
## ternaryplot {vcd}
require(vcd)
require(Zelig)
ternaryplot(x,
grid=FALSE,
dimnames.position = none,
pch=pch, col=col,
scale=1, main=NULL,
prop.size=FALSE
)
oldpar - par(no.readonly=TRUE)
ticlength - 0.01
## now the bottom internal ticks
x1-seq(0.1,0.9,by=0.1)
x2-x1
y1-rep(0,9)
y2-rep(ticlength,9)
segments(x1,y1,x2,y2)
text(x1,y1-0.03,as.character(rev(seq(10,90,by=10 #, cex=0.8)
## now the left internal ticks
y1-x1*sin60
x1-x1*0.5
x2-x1+ticlength*sin60
y2-y1-ticlength*0.5
segments(x1,y1,x2,y2)
text(x1-0.03,y1+0.015,as.character(seq(10,90,by=10)))
## now the right internal ticks
x1-rev(x1+0.5-ticlength*sin60)
x2-x1+ticlength*sin60
segments(x1,y2,x2,y1)
text(x2+0.03,y1+0.015,as.character(rev(seq(10,90,by=10
## the labels at the corners
par(xpd=TRUE)
# text(0.5,0.9,100% clay)
# text(-0.1,0,100% sand)
# text(1.1,0,100% loam)
## the axis labels
text(0.09,0.43,% Clay)
text(0.90,0.43,% Silt)
text(0.5,-0.1,% Sand)
# boundary of clay with extensions
x1-c(0.275,0.355,0.6)
x2-c(0.415,0.8,0.7)
y1-c(0.55*sin60,0.4*sin60,0.4*sin60)
y2-c(0.285*sin60,0.4*sin60,0.6*sin60)
segments(x1,y1,x2,y2, col=grey)
text(0.5,0.57,Clay, col=grey)
# lower bound of clay loam silty divider
x1-c(0.415,0.66)
x2-c(0.856,0.6)
y1-c(0.285*sin60,0.285*sin60)
y2-c(0.285*sin60,0.40*sin60)
segments(x1,y1,x2,y2, col=grey)
text(0.7,0.49*sin60,Silty, col=grey)
text(0.7,0.44*sin60,clay, col=grey)
text(0.72,0.36*sin60,Silty clay, col=grey)
text(0.73,0.32*sin60,loam, col=grey)
text(0.5,0.35*sin60,Clay loam, col=grey)
x1-c(0.185,0.1,0.37)
x2-c(0.37,0.37,0.415)
y1-c(0.37*sin60,0.2*sin60,0.2*sin60)
y2-c(0.37*sin60,0.2*sin60,0.285*sin60)
segments(x1,y1,x2,y2, col=grey)
text(0.28,0.43*sin60,Sandy, col=grey)
text(0.27,0.39*sin60,clay, col=grey)
text(0.27,0.3*sin60,Sandy clay, col=grey)
text(0.27,0.26*sin60,loam, col=grey)
# sand corner
x1-c(0.05,0.075)
x2-c(0.15,0.3)
y1-c(0.1*sin60,0.15*sin60)
y2-c(0,0)
segments(x1,y1,x2,y2, col=grey)
text(0.25,0.13*sin60,Sandy loam, col=grey)
text(0.14,0.07*sin60,Loamy, col=grey)
text(0.18,0.03*sin60,sand, col=grey)
text(0.06,0.021,Sand, col=grey)
x1-c(0.37,0.435,0.5,0.8,0.86)
x2-c(0.435,0.537,0.64,0.86,0.94)
y1-c(0.2*sin60,0.08*sin60,0,0,0.12*sin60)
y2-c(0.08*sin60,0.08*sin60,0.285*sin60,0.12*sin60,0.12*sin60)
segments(x1,y1,x2,y2, col=grey)
text(0.49,0.18*sin60,Loam, col=grey)
text(0.72,0.15*sin60,Silt loam, col=grey)
text(0.9,0.06*sin60,Silt, col=grey)
ternarypoints(x, pch = pch, col = col)
par(oldpar)
}
tmp - array(dim=c(10,3))
tmp[,2] - abs(rnorm(10)*20)
tmp[,3] - abs(rnorm(10)*10)
tmp[,1] - 100-tmp[,2]-tmp[,3]
col - rep(black,10)
pch - rep(1, 10)
plot.soiltexture(tmp,pch,col=black)
Sander Oom wrote:
Right,
Got the data points plotted on top of the soil texture background,
thanks to Jim and ternaryplot{vcd}! See code below.
Now there is some fine tuning to do, as it should really look like this
graph:
http://soil.scijournals.org/content/vol65/issue4/images/large/1038f2.jpeg
Things to do:
- rotate axis labels;
- correct small errors in class divisions;
- correct the partial covering of the bottom tick labels;
- rotate ticks in order to simplify viewing the graph.
Any help still appreciated!
Cheers,
Sander.
soil.triangle - function() {
oldpar - par(no.readonly=TRUE)
## now the bottom internal ticks
x1-seq(0.1,0.9,by=0.1)
x2-x1
y1-rep(0,9)
y2-rep(-0.02,9)
segments(x1,y1,x2,y2)
text(x1,y1-0.03,as.character(rev(seq(10,90,by=10 #, cex=0.8)
## now the left internal ticks
y1-x1*sin60
x1-x1*0.5
x2-x1+0.02*sin60
y2-y1-0.02*0.5
segments(x1,y1,x2,y2)
text(x1-0.03,y1+0.015,as.character(seq(10,90,by=10)))
## now the right internal ticks
x1-rev(x1+0.5-0.02*sin60)
x2-x1+0.02*sin60
segments(x1,y2,x2,y1)
text(x2+0.03,y1+0.015,as.character(rev(seq(10,90,by=10
## the labels at the corners
par(xpd=TRUE)
# text(0.5,0.9,100% clay)
# text(-0.1,0,100% sand)
# text(1.1,0,100% loam)
text(0.09,0.43,% Clay)
text(0.90,0.43,% Silt)
text(0.5,-0.1,% Sand)
# boundary of clay with extensions
x1-c(0.275,0.35,0.6)
x2-c(0.4,0.79,0.7)
y1-c(0.55*sin60,0.41*sin60,0.41*sin60)
y2-c(0.285*sin60,0.41*sin60,0.6*sin60)
segments(x1,y1,x2,y2, col=grey)
text(0.5,0.57,Clay, col=grey)
# lower bound of clay loam silty divider
x1-c
[R] How to get special (Hershey) font symbols into plot axis labels? [Revisited]
Dear R users, I would like to use sub- and super-script in axis labels. I assume this is best done using Hershey symbols. When trying to find information on using Hershey font symbols in axis labels, I came across the following discussion thread: http://maths.newcastle.edu.au/~rking/R/help/02a/1857.html Have Hershey font implementations moved on since then? Thanks, Sander. -- Dr Sander P. Oom Animal, Plant and Environmental Sciences, University of the Witwatersrand Private Bag 3, Wits 2050, South Africa Tel (work) +27 (0)11 717 64 04 Tel (home) +27 (0)18 297 44 51 Fax +27 (0)18 299 24 64 Email [EMAIL PROTECTED] Web www.oomvanlieshout.net/sander __ 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
Re: [R] Plot range resizing when adding additiona lines
Hikel, Jerry wrote:
Hi -- I have searched the documentation and archives on graphing
capabilities in R for the past couple of hours, but I haven't been able
to find anything directly related to my problem.
I want to create a plot with several lines displayed on it. I want each
line to be displayed in a different color, so it seems the only way to
do this is to create a plot and then use the lines() function. However,
once i create the original plot, when I add new lines to the plot, the
axes range stays the same, and the new lines often extend beyond the
range of the original axes, cutting off the data of the added lines.
Is there any way to force the plot to resize it's axes range when new
lines are added, or is there some other optional way of implementing
this? Thanks.
DISCLAIMER: This e-mail message and any attachments are inte...{{dropped}}
__
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
You can check the range of all items to plot and take the minimum and
maximum value across all items. Then draw the plot with the axis range
explicitely set to 'minimum'-'maximum'!
That is assuming you know the items you are drawing before hand.
Sander.
--
Dr Sander P. Oom
Animal, Plant and Environmental Sciences,
University of the Witwatersrand
Private Bag 3, Wits 2050, South Africa
Tel (work) +27 (0)11 717 64 04
Tel (home) +27 (0)18 297 44 51
Fax +27 (0)18 299 24 64
Email [EMAIL PROTECTED]
Web www.oomvanlieshout.net/sander
__
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
Re: [R] adjusted p-values with TukeyHSD?
Hi Chris and Chris, I was keeping my eye on this thread as I have also been discovering multiple comparisons recently. Your instructions are very clear! Thanks. Now I would love to see an R boffin write a nifty function to produce a graphical representation of the multiple comparison, like this one: http://www.theses.ulaval.ca/2003/21026/21026024.jpg Should not be too difficult.[any one up for the challenge?] I came across more multiple comparison info here; http://www.agr.kuleuven.ac.be/vakken/statisticsbyR/ANOVAbyRr/multiplecomp.htm Cheers, Sander. Christoph Buser wrote: Dear Christoph You can use the multcomp package. Please have a look at the following example: library(multcomp) The first two lines were already proposed by Erin Hodgess: summary(fm1 - aov(breaks ~ wool + tension, data = warpbreaks)) TukeyHSD(fm1, tension, ordered = TRUE) Tukey multiple comparisons of means 95% family-wise confidence level factor levels have been ordered Fit: aov(formula = breaks ~ wool + tension, data = warpbreaks) $tension difflwr upr M-H 4.72 -4.6311985 14.07564 L-H 14.72 5.3688015 24.07564 L-M 10.00 0.6465793 19.35342 By using the functions simtest or simint you can get the p-values, too: summary(simtest(breaks ~ wool + tension, data = warpbreaks, whichf=tension, type = Tukey)) Simultaneous tests: Tukey contrasts Call: simtest.formula(formula = breaks ~ wool + tension, data = warpbreaks, whichf = tension, type = Tukey) Tukey contrasts for factor tension, covariable: wool Contrast matrix: tensionL tensionM tensionH tensionM-tensionL 0 0 -110 tensionH-tensionL 0 0 -101 tensionH-tensionM 0 00 -11 Absolute Error Tolerance: 0.001 Coefficients: Estimate t value Std.Err. p raw p Bonf p adj tensionH-tensionL -14.722 -3.8023.872 0.000 0.001 0.001 tensionM-tensionL -10.000 -2.5823.872 0.013 0.026 0.024 tensionH-tensionM -4.722 -1.2193.872 0.228 0.228 0.228 or if you prefer to get the confidence intervals, too, you can use: summary(simint(breaks ~ wool + tension, data = warpbreaks, whichf=tension, type = Tukey)) Simultaneous 95% confidence intervals: Tukey contrasts Call: simint.formula(formula = breaks ~ wool + tension, data = warpbreaks, whichf = tension, type = Tukey) Tukey contrasts for factor tension, covariable: wool Contrast matrix: tensionL tensionM tensionH tensionM-tensionL 0 0 -110 tensionH-tensionL 0 0 -101 tensionH-tensionM 0 00 -11 Absolute Error Tolerance: 0.001 95 % quantile: 2.415 Coefficients: Estimate 2.5 % 97.5 % t value Std.Err. p raw p Bonf p adj tensionM-tensionL -10.000 -19.352 -0.648 -2.5823.872 0.013 0.038 0.034 tensionH-tensionL -14.722 -24.074 -5.370 -3.8023.872 0.000 0.001 0.001 tensionH-tensionM -4.722 -14.074 4.630 -1.2193.872 0.228 0.685 0.447 - Please be careful: The resulting confidence intervals in simint are not associated with the p-values from 'simtest' as it is described in the help page of the two functions. - I had not the time to check the differences in the function or read the references given on the help page. If you are interested in the function you can check those to find out which one you prefer. Best regards, Christoph Buser -- Christoph Buser [EMAIL PROTECTED] Seminar fuer Statistik, LEO C13 ETH (Federal Inst. Technology) 8092 Zurich SWITZERLAND phone: x-41-44-632-4673 fax: 632-1228 http://stat.ethz.ch/~buser/ -- Christoph Strehblow writes: hi list, i have to ask you again, having tried and searched for several days... i want to do a TukeyHSD after an Anova, and want to get the adjusted p-values after the Tukey Correction. i found the p.adjust function, but it can only correct for holm, hochberg, bonferroni, but not Tukey. Is it not possbile to get adjusted p-values after Tukey-correction? sorry, if this is an often-answered-question, but i didn´t find it on the list archive... thx a lot, list, Chris Christoph Strehblow, MD Department of Rheumatology, Diabetes and Endocrinology Wilhelminenspital, Vienna, Austria [EMAIL PROTECTED] __ 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 __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read
Re: [R] adjusted p-values with TukeyHSD?
Shame I can not get hold of Hsu, J. C. and M. Peruggia (1994) just now.
I am quite curious to see what their graphs look like. Would you be able
to give an example in R.? ;-)
The graph I put forward is typically used by ecologists to summarize
data. It comes down to a simple means plot with error bars. Significant
differences of multiple comparisons are then added using the letters a,
b, c etc. If two bars have the same letter, they are not significantly
different. It can become quite complicated when mean one is different
from mean three but not from mean two and mean two is different from
mean three but not mean one. You then get: a, ab, c for mean one, two
and three respectively.
Of course what is often used does not constitute the best way of doing it.
Sander.
Liaw, Andy wrote:
From: Sander Oom
Hi Chris and Chris,
I was keeping my eye on this thread as I have also been discovering
multiple comparisons recently. Your instructions are very
clear! Thanks.
One thing to note, though: Multcomp does not do Dunnett's or
Tukey's multiple comparisons per se. Those names in multcomp
refer to the contrasts being used (comparison to a control for
Dunnett and all pairwise comparison for Tukey). The actual
methods used are as described in the references of the help
pages.
Now I would love to see an R boffin write a nifty function to
produce a
graphical representation of the multiple comparison, like this one:
http://www.theses.ulaval.ca/2003/21026/21026024.jpg
Should not be too difficult.[any one up for the challenge?]
I beg to differ: That's probably as bad a way as one can use to
graphically show multiple comparison. The shaded bars serve no
purpose.
Two alternatives that I'm aware of are
- Multiple comparison circles, due to John Sall, and not
surprisingly, implemented in JMP and SAS/Insight. See:
http://support.sas.com/documentation/onlinedoc/v7/whatsnew/insight/sect4.htm
- The mean-mean display proposed by Hsu and Peruggia:
Hsu, J. C. and M. Peruggia (1994).
Graphical representations of Tukey's multiple comparison method.
Journal of Computational and Graphical Statistics 3, 143{161
Andy
I came across more multiple comparison info here;
http://www.agr.kuleuven.ac.be/vakken/statisticsbyR/ANOVAbyRr/m
ultiplecomp.htm
Cheers,
Sander.
Christoph Buser wrote:
Dear Christoph
You can use the multcomp package. Please have a look at the
following example:
library(multcomp)
The first two lines were already proposed by Erin Hodgess:
summary(fm1 - aov(breaks ~ wool + tension, data = warpbreaks))
TukeyHSD(fm1, tension, ordered = TRUE)
Tukey multiple comparisons of means
95% family-wise confidence level
factor levels have been ordered
Fit: aov(formula = breaks ~ wool + tension, data = warpbreaks)
$tension
difflwr upr
M-H 4.72 -4.6311985 14.07564
L-H 14.72 5.3688015 24.07564
L-M 10.00 0.6465793 19.35342
By using the functions simtest or simint you can get the
p-values, too:
summary(simtest(breaks ~ wool + tension, data = warpbreaks,
whichf=tension,
type = Tukey))
Simultaneous tests: Tukey contrasts
Call:
simtest.formula(formula = breaks ~ wool + tension, data =
warpbreaks,
whichf = tension, type = Tukey)
Tukey contrasts for factor tension, covariable: wool
Contrast matrix:
tensionL tensionM tensionH
tensionM-tensionL 0 0 -110
tensionH-tensionL 0 0 -101
tensionH-tensionM 0 00 -11
Absolute Error Tolerance: 0.001
Coefficients:
Estimate t value Std.Err. p raw p Bonf p adj
tensionH-tensionL -14.722 -3.8023.872 0.000 0.001 0.001
tensionM-tensionL -10.000 -2.5823.872 0.013 0.026 0.024
tensionH-tensionM -4.722 -1.2193.872 0.228 0.228 0.228
or if you prefer to get the confidence intervals, too, you can
use:
summary(simint(breaks ~ wool + tension, data = warpbreaks,
whichf=tension,
type = Tukey))
Simultaneous 95% confidence intervals: Tukey contrasts
Call:
simint.formula(formula = breaks ~ wool + tension, data =
warpbreaks,
whichf = tension, type = Tukey)
Tukey contrasts for factor tension, covariable: wool
Contrast matrix:
tensionL tensionM tensionH
tensionM-tensionL 0 0 -110
tensionH-tensionL 0 0 -101
tensionH-tensionM 0 00 -11
Absolute Error Tolerance: 0.001
95 % quantile: 2.415
Coefficients:
Estimate 2.5 % 97.5 % t value Std.Err.
p raw p Bonf p adj
tensionM-tensionL -10.000 -19.352 -0.648 -2.5823.872
0.013 0.038 0.034
tensionH-tensionL -14.722 -24.074 -5.370 -3.8023.872
0.000 0.001 0.001
tensionH-tensionM -4.722 -14.074 4.630 -1.2193.872
0.228 0.685 0.447
-
Please be careful: The resulting confidence intervals in
simint are not associated with the p-values
Re: [R] Conflict between xtable and Hmisc when using Sweave?
Dear David,
I would like to use summarize(Hmisc) and print.xtable(xtable) in a
single Sweave document, but a conflict with the 'label' function
prohibits this at the moment!
Would you be able to correct the conflicting code? I will gladly test
the new package!
I have tried latex(Hmisc) to export the anova table, but results are not
promising! I prefer xtable!!
Thanks,
Sander.
Frank E Harrell Jr wrote:
Sander Oom wrote:
Dear Frank,
I have a Sweave document in which I export anova (aov) tables to Latex
and calculate some summary statistics with summarize{Hmisc} for a
graph (as in the example below).
I currently use the following code for the aov tables:
results=tex=
tmp - datGrassHC[datGrassHC$Loc 0 datGrassHC$Loc 9 ,]
tmpAov - aov(Height~Geology*Altitude*Origin*BinInOut , data=tmp)
tmpTable - xtable (tmpAov ,
caption=ANOVA table for vegetation height.,
label=tab:AnovaHeight
)
print.xtable(tmpTable, type=latex, floating=TRUE,
table.placement=ht, caption.placement=top,
latex.environments=c(center))
)
@
I used xtables, because it has a working aov example. I would be happy
to use an alternative if I knew how! Would you have sample code to
illustrate how to export an aov table to Latex using latex{Hmisc}.
Thanks very much for your help,
Sander.
Frank E Harrell Jr wrote:
Sander Oom wrote:
Dear R users,
The Sweave code below runs fine, as it is. However, an error occurs
when the line 'library(xtable)' is uncommented:
Error: chunk 1
Error in label-(`*tmp*`, value = month) :
no applicable method for label-
Is anybody aware of this and knows a workaround?
Thanks,
Sander.
***
\documentclass[a4paper]{article}
\title{Sweave Test for summarize}
\author{Sander Oom}
\usepackage{a4wide}
\begin{document}
\maketitle
\begin{figure}[ht]
\begin{center}
fig=TRUE,echo=FALSE=
# library(xtable)
library(Hmisc)
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)
print(xYplot(Cbind(y,Lower,Upper) ~ month, groups=year, data=s,
keys='lines', method='alt', type='b'))
@
\end{center}
\end{figure}
\end{document}
version
_
platform i686-pc-linux-gnu
arch i686
os linux-gnu
system i686, linux-gnu
status
major2
minor1.0
year 2005
month04
day 18
language R
I feel this is an xtable problem because Hmisc has being using label
and label- since 1991.
Frank
There are ways to make functions from one area override those from
another, but the real solution is to ask the xtable author not to have
functions that conflict with the (older) Hmisc package. -Frank
--
Dr Sander P. Oom
Animal, Plant and Environmental Sciences,
University of the Witwatersrand
Private Bag 3, Wits 2050, South Africa
Tel (work) +27 (0)11 717 64 04
Tel (home) +27 (0)18 297 44 51
Fax +27 (0)18 299 24 64
Email [EMAIL PROTECTED]
Web www.oomvanlieshout.net/sander
__
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
Re: [R] Conflict between xtable and Hmisc when using Sweave?
Hi Andy and Gabor,
Thanks for your help so far! I am discovering another R dimension.
Trying to put my head around all thisthe conflict actually exposes
itself when calling summarize(Hmisc). Summarize(Hmisc) calls label
internally, so I can not call it explicitly. Simply calling
label(xtable) explicitly will not solve the problem with summarize(Hmisc).
Thus, I should use namespaces as Andy is suggesting. Now I just need to
know how I 'add namespace' to a library? Does 'loadNamespace' have
something to do with it?
Thanks very much for your help!
Sander.
## From Venables and Ripley (2002) p.165.
N - c(0,1,0,1,1,1,0,0,0,1,1,0,1,1,0,0,1,0,1,0,1,1,0,0)
P - c(1,1,0,0,0,1,0,1,1,1,0,0,0,1,0,1,1,0,0,1,0,1,1,0)
K - c(1,0,0,1,0,1,1,0,0,1,0,1,0,1,1,0,0,0,1,1,1,0,1,0)
yield -c(49.5,62.8,46.8,57.0,59.8,58.5,55.5,56.0,
62.8,55.8,69.5,55.0,
62.0,48.8,45.5,44.2,52.0,
51.5,49.8,48.8,57.2,59.0,53.2,56.0)
npk - data.frame(block=gl(6,4), N=factor(N), P=factor(P),
K=factor(K), yield=yield)
## to show the effects of re-ordering terms contrast the two fits
tmpAov - aov(yield ~ block + N * P + K, npk)
tmpTable - xtable(tmpAov , caption=Test export of ANOVA table.,
label=tab:Anova)
print.xtable(tmpTable, type=latex, floating=TRUE,
table.placement=h, caption.placement=top,
latex.environments=c(center))
Alternatively, using namespace for xtable:
tmpTable - xtable(tmpAov , caption=Test export of ANOVA table.)
xtable:::label(tmpTable) - paste(tab:Anova)
print.xtable(tmpTable, type=latex, floating=TRUE,
table.placement=ht, caption.placement=top,
latex.environments=c(center))
Gabor Grothendieck wrote:
Even without a namespace one could explicitly reference the label
in xtable via:
xtable.label - get(label, package:xtable)
On 5/16/05, Liaw, Andy [EMAIL PROTECTED] wrote:
One possible solution without renaming the functions is to add namespace to
either xtable or Hmisc. Given the size of Hmisc, it probably would be much
easier to do that with xtable.
With namespace in xtable, you can do xtable:::label() to refer to the
label() in xtable specifically.
Andy
From: Of Sander Oom
Dear David,
I would like to use summarize(Hmisc) and print.xtable(xtable) in a
single Sweave document, but a conflict with the 'label' function
prohibits this at the moment!
Would you be able to correct the conflicting code? I will gladly test
the new package!
I have tried latex(Hmisc) to export the anova table, but
results are not
promising! I prefer xtable!!
Thanks,
Sander.
Frank E Harrell Jr wrote:
Sander Oom wrote:
Dear Frank,
I have a Sweave document in which I export anova (aov)
tables to Latex
and calculate some summary statistics with summarize{Hmisc} for a
graph (as in the example below).
I currently use the following code for the aov tables:
results=tex=
tmp - datGrassHC[datGrassHC$Loc 0 datGrassHC$Loc 9 ,]
tmpAov - aov(Height~Geology*Altitude*Origin*BinInOut , data=tmp)
tmpTable - xtable (tmpAov ,
caption=ANOVA table for vegetation height.,
label=tab:AnovaHeight
)
print.xtable(tmpTable, type=latex, floating=TRUE,
table.placement=ht, caption.placement=top,
latex.environments=c(center))
)
@
I used xtables, because it has a working aov example. I
would be happy
to use an alternative if I knew how! Would you have sample code to
illustrate how to export an aov table to Latex using latex{Hmisc}.
Thanks very much for your help,
Sander.
Frank E Harrell Jr wrote:
Sander Oom wrote:
Dear R users,
The Sweave code below runs fine, as it is. However, an
error occurs
when the line 'library(xtable)' is uncommented:
Error: chunk 1
Error in label-(`*tmp*`, value = month) :
no applicable method for label-
Is anybody aware of this and knows a workaround?
Thanks,
Sander.
***
\documentclass[a4paper]{article}
\title{Sweave Test for summarize}
\author{Sander Oom}
\usepackage{a4wide}
\begin{document}
\maketitle
\begin{figure}[ht]
\begin{center}
fig=TRUE,echo=FALSE=
# library(xtable)
library(Hmisc)
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)
print(xYplot(Cbind(y,Lower,Upper) ~ month, groups=year, data=s,
keys='lines', method='alt', type='b'))
@
\end{center}
\end{figure}
\end{document}
version
_
platform i686-pc-linux-gnu
arch i686
os linux-gnu
system i686, linux-gnu
status
major2
minor1.0
year 2005
month04
day 18
language R
I feel this is an xtable problem because Hmisc has being
using label
and label- since 1991.
Frank
There are ways to make functions from one area override those from
another, but the real solution is to ask the xtable author
not to have
functions that conflict with the (older) Hmisc package. -Frank
--
Dr Sander P. Oom
Animal, Plant and Environmental
Re: [R] Conflict between xtable and Hmisc when using Sweave?
I tried to follow your suggestions, but without success:
Error: couldn't find function xtable.mylabel-
... resulting from the code below.
Any suggestions?
Thanks,
Sander.
library(xtable)
xtable.mylabel - get(label, package:xtable)
library(Hmisc) # provides summarize
set.seed(1)
temperature - rnorm(300, 70, 10)
month - sample(1:12, 300, TRUE)
year - sample(2000:2001, 300, TRUE)
g - function(x)c(Mean=mean(x,na.rm=TRUE),Median=median(x,na.rm=TRUE))
summarize(temperature, month, g)
## From Venables and Ripley (2002) p.165.
N - c(0,1,0,1,1,1,0,0,0,1,1,0,1,1,0,0,1,0,1,0,1,1,0,0)
P - c(1,1,0,0,0,1,0,1,1,1,0,0,0,1,0,1,1,0,0,1,0,1,1,0)
K - c(1,0,0,1,0,1,1,0,0,1,0,1,0,1,1,0,0,0,1,1,1,0,1,0)
yield -c(49.5,62.8,46.8,57.0,59.8,58.5,55.5,56.0,
62.8,55.8,69.5,55.0,
62.0,48.8,45.5,44.2,52.0,
51.5,49.8,48.8,57.2,59.0,53.2,56.0)
npk - data.frame(block=gl(6,4), N=factor(N), P=factor(P),
K=factor(K), yield=yield)
## to show the effects of re-ordering terms contrast the two fits
tmpAov - aov(yield ~ block + N * P + K, npk)
tmpTable - xtable (tmpAov ,
caption=ANOVA table for vegetation height.)
xtable.mylabel(tmpTable) - paste(tab:AnovaHeight)
print.xtable(tmpTable, type=latex, floating=TRUE,
table.placement=h, caption.placement=top,
latex.environments=c(center),
title=first.word(deparse(substitute(object))),
append=FALSE
)
Liaw, Andy wrote:
You need to add the namespace to the source package, by adding a NAMESPACE
file. There's an R News article by Prof. Tierney on how to do this. Also
see the `Writing R Extensions' manual. You should get the package
maintainer to do that, as that constitute a change in the package source
code.
Short of that, you should make sure that Hmisc is loaded later than xtable,
and use something like what Gabor suggested to access label() in xtable. (I
would use some other name, though: label() in xtable is already an S3
generic).
Andy
From: Sander Oom
Hi Andy and Gabor,
Thanks for your help so far! I am discovering another R dimension.
Trying to put my head around all thisthe conflict
actually exposes
itself when calling summarize(Hmisc). Summarize(Hmisc) calls label
internally, so I can not call it explicitly. Simply calling
label(xtable) explicitly will not solve the problem with
summarize(Hmisc).
Thus, I should use namespaces as Andy is suggesting. Now I
just need to
know how I 'add namespace' to a library? Does 'loadNamespace' have
something to do with it?
Thanks very much for your help!
Sander.
## From Venables and Ripley (2002) p.165.
N - c(0,1,0,1,1,1,0,0,0,1,1,0,1,1,0,0,1,0,1,0,1,1,0,0)
P - c(1,1,0,0,0,1,0,1,1,1,0,0,0,1,0,1,1,0,0,1,0,1,1,0)
K - c(1,0,0,1,0,1,1,0,0,1,0,1,0,1,1,0,0,0,1,1,1,0,1,0)
yield -c(49.5,62.8,46.8,57.0,59.8,58.5,55.5,56.0,
62.8,55.8,69.5,55.0,
62.0,48.8,45.5,44.2,52.0,
51.5,49.8,48.8,57.2,59.0,53.2,56.0)
npk - data.frame(block=gl(6,4), N=factor(N), P=factor(P),
K=factor(K), yield=yield)
## to show the effects of re-ordering terms contrast the two fits
tmpAov - aov(yield ~ block + N * P + K, npk)
tmpTable - xtable(tmpAov , caption=Test export of ANOVA table.,
label=tab:Anova)
print.xtable(tmpTable, type=latex, floating=TRUE,
table.placement=h, caption.placement=top,
latex.environments=c(center))
Alternatively, using namespace for xtable:
tmpTable - xtable(tmpAov , caption=Test export of ANOVA table.)
xtable:::label(tmpTable) - paste(tab:Anova)
print.xtable(tmpTable, type=latex, floating=TRUE,
table.placement=ht, caption.placement=top,
latex.environments=c(center))
Gabor Grothendieck wrote:
Even without a namespace one could explicitly reference the label
in xtable via:
xtable.label - get(label, package:xtable)
On 5/16/05, Liaw, Andy [EMAIL PROTECTED] wrote:
One possible solution without renaming the functions is to
add namespace to
either xtable or Hmisc. Given the size of Hmisc, it
probably would be much
easier to do that with xtable.
With namespace in xtable, you can do xtable:::label() to
refer to the
label() in xtable specifically.
Andy
From: Of Sander Oom
Dear David,
I would like to use summarize(Hmisc) and print.xtable(xtable) in a
single Sweave document, but a conflict with the 'label' function
prohibits this at the moment!
Would you be able to correct the conflicting code? I will
gladly test
the new package!
I have tried latex(Hmisc) to export the anova table, but
results are not
promising! I prefer xtable!!
Thanks,
Sander.
Frank E Harrell Jr wrote:
Sander Oom wrote:
Dear Frank,
I have a Sweave document in which I export anova (aov)
tables to Latex
and calculate some summary statistics with summarize{Hmisc} for a
graph (as in the example below).
I currently use the following code for the aov tables:
results=tex=
tmp - datGrassHC[datGrassHC$Loc 0 datGrassHC$Loc 9 ,]
tmpAov - aov(Height~Geology*Altitude*Origin*BinInOut
, data=tmp)
tmpTable - xtable (tmpAov ,
caption=ANOVA table for vegetation height.,
label=tab:AnovaHeight
)
print.xtable
Re: [R] Conflict between xtable and Hmisc when using Sweave?
Hi David,
Thanks for creating and supporting xtable. Glad I can contribute by
pointing out problems!
An earlier response from Frank Harrell Jr. suggests that the
responsibility to resolve the conflict lies primarily with you, as Hmisc
'was there first'.
Not sure how disputes over package conflicts are generally resolved!?
I'll await a final decision!
Thanks again,
Sander.
David B. Dahl wrote:
Sander,
Thanks for pointing out the conflict between Hmisc and xtable. I am not
sure I have a good solution.
My understand of the namespace solution is that packages can specify
which variables to export for use by the package users. The label
function is not an internal function, rather one what is intended for
the user.
Renaming the label() would resolve the conflict with the Hmisc package,
but make xtable not compatible with previous versions.
As noted by Andy, label() in xtable is an S3 generic with an
implementation label.xtable() specific to the xtable package. Perhaps
Frank of Hmisc might be willing to make his follow the S3 generic naming
convention?
I am open to suggestions and, more especially, code. Thanks for using
xtable.
-- David
Liaw, Andy wrote:
You need to add the namespace to the source package, by adding a
NAMESPACE
file. There's an R News article by Prof. Tierney on how to do this.
Also
see the `Writing R Extensions' manual. You should get the package
maintainer to do that, as that constitute a change in the package source
code.
Short of that, you should make sure that Hmisc is loaded later than
xtable,
and use something like what Gabor suggested to access label() in
xtable. (I
would use some other name, though: label() in xtable is already an S3
generic).
Andy
From: Sander Oom
Hi Andy and Gabor,
Thanks for your help so far! I am discovering another R dimension.
Trying to put my head around all thisthe conflict actually
exposes itself when calling summarize(Hmisc). Summarize(Hmisc) calls
label internally, so I can not call it explicitly. Simply calling
label(xtable) explicitly will not solve the problem with
summarize(Hmisc).
Thus, I should use namespaces as Andy is suggesting. Now I just need
to know how I 'add namespace' to a library? Does 'loadNamespace' have
something to do with it?
Thanks very much for your help!
Sander.
## From Venables and Ripley (2002) p.165.
N - c(0,1,0,1,1,1,0,0,0,1,1,0,1,1,0,0,1,0,1,0,1,1,0,0)
P - c(1,1,0,0,0,1,0,1,1,1,0,0,0,1,0,1,1,0,0,1,0,1,1,0)
K - c(1,0,0,1,0,1,1,0,0,1,0,1,0,1,1,0,0,0,1,1,1,0,1,0)
yield -c(49.5,62.8,46.8,57.0,59.8,58.5,55.5,56.0,
62.8,55.8,69.5,55.0,
62.0,48.8,45.5,44.2,52.0,
51.5,49.8,48.8,57.2,59.0,53.2,56.0)
npk - data.frame(block=gl(6,4), N=factor(N), P=factor(P),
K=factor(K), yield=yield)
## to show the effects of re-ordering terms contrast the two fits
tmpAov - aov(yield ~ block + N * P + K, npk)
tmpTable - xtable(tmpAov , caption=Test export of ANOVA table.,
label=tab:Anova)
print.xtable(tmpTable, type=latex, floating=TRUE,
table.placement=h, caption.placement=top,
latex.environments=c(center))
Alternatively, using namespace for xtable:
tmpTable - xtable(tmpAov , caption=Test export of ANOVA table.)
xtable:::label(tmpTable) - paste(tab:Anova)
print.xtable(tmpTable, type=latex, floating=TRUE,
table.placement=ht, caption.placement=top,
latex.environments=c(center))
Gabor Grothendieck wrote:
Even without a namespace one could explicitly reference the label
in xtable via:
xtable.label - get(label, package:xtable)
On 5/16/05, Liaw, Andy [EMAIL PROTECTED] wrote:
One possible solution without renaming the functions is to
add namespace to
either xtable or Hmisc. Given the size of Hmisc, it
probably would be much
easier to do that with xtable.
With namespace in xtable, you can do xtable:::label() to
refer to the
label() in xtable specifically.
Andy
From: Of Sander Oom
Dear David,
I would like to use summarize(Hmisc) and print.xtable(xtable) in a
single Sweave document, but a conflict with the 'label' function
prohibits this at the moment!
Would you be able to correct the conflicting code? I will
gladly test
the new package!
I have tried latex(Hmisc) to export the anova table, but
results are not
promising! I prefer xtable!!
Thanks,
Sander.
Frank E Harrell Jr wrote:
Sander Oom wrote:
Dear Frank,
I have a Sweave document in which I export anova (aov)
tables to Latex
and calculate some summary statistics with summarize{Hmisc} for a
graph (as in the example below).
I currently use the following code for the aov tables:
results=tex=
tmp - datGrassHC[datGrassHC$Loc 0 datGrassHC$Loc 9 ,]
tmpAov - aov(Height~Geology*Altitude*Origin*BinInOut
, data=tmp)
tmpTable - xtable (tmpAov ,
caption=ANOVA table for vegetation height.,
label=tab:AnovaHeight
)
print.xtable(tmpTable, type=latex, floating=TRUE,
table.placement=ht, caption.placement=top
Re: [R] Problem with data frame when using xYplot?
Problem solved! I was so focused on reproducing the plotmeans() functionality with xYplot() that I completely overlooked the fact that my data does not allow a x-y plot, as only Sodium is a numeric variable while Position and AltGeo are factors! Using unclass() to make Position a numeric variable does the trick: tmp$Position - unclass(tmp$Position) The code below does the trick. Now I only need to figure out how to tweak the x axis to pretend I am plotting a factor, i.e. plotting labels Inside and Outside. Cheers, Sander. library(Hmisc) library(Lattice) tmp - structure(list(Position = structure(as.integer(c(1, 2, 1, 2, 1, 2, 1, 2)), .Label = c(Inside, Outside), class = factor), AltGeo = structure(as.integer(c(1, 1, 2, 2, 3, 3, 4, 4)), .Label = c(Basalt-High, Basalt-Low, Quartz-High, Quartz-Low), class = factor), Sodium = c(27.3, 26.9, 25, 18.1, 4.67, 5.56, 10.7, 5.67 ), SD = c(5.3851648071345, 2.42097317438899, 20.1618451536560, 15.2679766541317, 5.45435605731786, 8.09492296305393, 10.6183802907976, 8.06225774829855), Nobs = c(9, 9, 9, 9, 9, 9, 9, 9), Lower = c(25.5382783976218, 26.0818978307592, 18.2793849487813, 13.0217855597339, 2.84854798089405, 2.85724790120425, 7.12720656973412, 2.97924741723382), Upper = c(29.1283882690448, 27.6958799470186, 31.7206150512187, 23.2004366624884, 6.48478535243929, 8.25386320990686, 14.2061267635992, 8.35408591609952)), .Names = c(Position, AltGeo, Sodium, SD, Nobs, Lower, Upper), row.names = c(1, 2, 3, 4, 5, 6, 7, 8), class = data.frame) tmp$Position - unclass(tmp$Position) xYplot(Cbind(Sodium,Lower,Upper) ~ Position|AltGeo, groups=AltGeo, data=tmp, ylim=c(min(tmp$Lower)-1,max(tmp$Upper)+1), xlab=Position, ylab=Sodium ) Sander Oom wrote: Dear all, I am trying to plot means and error bars using xYplot, but I get an error message from xYplot which I can not figure out: Error in Summary.factor(..., na.rm = na.rm) : range not meaningful for factors The data frame (tmpNa) was created using aggregate. I have used dump to created the code below, which generates the same error. Can anybody tell me what is wrong with the data frame? Thanks in advance, Sander. library(Hmisc) tmpNa - structure(list(Position = structure(as.integer(c(1, 2, 1, 2, 1, 2, 1, 2)), .Label = c(Inside, Outside), class = factor), AltGeo = structure(as.integer(c(1, 1, 2, 2, 3, 3, 4, 4)), .Label = c(Basalt-High, Basalt-Low, Quartz-High, Quartz-Low), class = factor), Sodium = c(27.3, 26.9, 25, 18.1, 4.67, 5.56, 10.7, 5.67 ), SD = c(5.3851648071345, 2.42097317438899, 20.1618451536560, 15.2679766541317, 5.45435605731786, 8.09492296305393, 10.6183802907976, 8.06225774829855), Nobs = c(9, 9, 9, 9, 9, 9, 9, 9), Lower = c(25.5382783976218, 26.0818978307592, 18.2793849487813, 13.0217855597339, 2.84854798089405, 2.85724790120425, 7.12720656973412, 2.97924741723382), Upper = c(29.1283882690448, 27.6958799470186, 31.7206150512187, 23.2004366624884, 6.48478535243929, 8.25386320990686, 14.2061267635992, 8.35408591609952)), .Names = c(Position, AltGeo, Sodium, SD, Nobs, Lower, Upper), row.names = c(1, 2, 3, 4, 5, 6, 7, 8), class = data.frame) xYplot(Cbind(Sodium,Lower,Upper) ~ AltGeo, groups=Position, data=tmpNa) version _ platform i686-pc-linux-gnu arch i686 os linux-gnu system i686, linux-gnu status major2 minor1.0 year 2005 month04 day 18 language R -- Dr. Sander P. Oom Animal, Plant and Environmental Sciences, University of the Witwatersrand Private Bag 3, Wits 2050, South Africa Tel (work) +27 (0)11 717 64 04 Tel (home) +27 (0)18 297 44 51 Fax +27 (0)18 299 24 64 Email [EMAIL PROTECTED] Web www.oomvanlieshout.net/sander __ 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
[R] Conflict between xtable and Hmisc when using Sweave?
Dear R users,
The Sweave code below runs fine, as it is. However, an error occurs when
the line 'library(xtable)' is uncommented:
Error: chunk 1
Error in label-(`*tmp*`, value = month) :
no applicable method for label-
Is anybody aware of this and knows a workaround?
Thanks,
Sander.
***
\documentclass[a4paper]{article}
\title{Sweave Test for summarize}
\author{Sander Oom}
\usepackage{a4wide}
\begin{document}
\maketitle
\begin{figure}[ht]
\begin{center}
fig=TRUE,echo=FALSE=
# library(xtable)
library(Hmisc)
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)
print(xYplot(Cbind(y,Lower,Upper) ~ month, groups=year, data=s,
keys='lines', method='alt', type='b'))
@
\end{center}
\end{figure}
\end{document}
version
_
platform i686-pc-linux-gnu
arch i686
os linux-gnu
system i686, linux-gnu
status
major2
minor1.0
year 2005
month04
day 18
language R
--
Dr. Sander P. Oom
Animal, Plant and Environmental Sciences,
University of the Witwatersrand
Private Bag 3, Wits 2050, South Africa
Tel (work) +27 (0)11 717 64 04
Tel (home) +27 (0)18 297 44 51
Fax +27 (0)18 299 24 64
Email [EMAIL PROTECTED]
Web www.oomvanlieshout.net/sander
__
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
Re: [R] Problem with data frame when using xYplot?
I have edited the code (hacked from another graph) to provide more
control over the different elements of the graph. Now we have a graph at
publication quality! Slowly the power of R graphics is shining through
the thick cloud of options! Beautiful.
library(Hmisc)
library(lattice)
ltheme - canonical.theme(color = FALSE) ## in-built BW theme
ltheme$strip.background$col - transparent ## change strip bg
lattice.options(default.theme = ltheme) ## set as default
tmp -
structure(list(Position = structure(as.integer(c(1, 2, 1, 2,
1, 2, 1, 2)), .Label = c(Inside, Outside), class = factor),
AltGeo = structure(as.integer(c(1, 1, 2, 2, 3, 3, 4, 4)), .Label =
c(Basalt-High,
Basalt-Low, Quartz-High, Quartz-Low), class = factor),
Sodium = c(27.3, 26.9, 25, 18.1,
4.67, 5.56, 10.7, 5.67
), SD = c(5.3851648071345, 2.42097317438899, 20.1618451536560,
15.2679766541317, 5.45435605731786, 8.09492296305393, 10.6183802907976,
8.06225774829855), Nobs = c(9, 9, 9, 9, 9, 9, 9, 9), Lower =
c(25.5382783976218,
26.0818978307592, 18.2793849487813, 13.0217855597339, 2.84854798089405,
2.85724790120425, 7.12720656973412, 2.97924741723382), Upper =
c(29.1283882690448,
27.6958799470186, 31.7206150512187, 23.2004366624884, 6.48478535243929,
8.25386320990686, 14.2061267635992, 8.35408591609952)), .Names =
c(Position,
AltGeo, Sodium, SD, Nobs, Lower, Upper), row.names = c(1,
2, 3, 4, 5, 6, 7, 8), class = data.frame)
tmp$xvar - rep(1:4, each=2)+rep(c(-.05,.05), 4)
tmp
sp - list(superpose.symbol = list(pch = c(16,1), cex = 1))
xYplot(Cbind(Sodium,Lower,Upper) ~ xvar, groups=Position, data=tmp,
scales=list(y='free',x=list(at=1:4, labels=levels(tmp$AltGeo))),
xlim=c(0.5, 4.5), ylim=c(min(tmp$Lower)-1,max(tmp$Upper)+1),
xlab='AltGeo', ylab='Sodium',
panel = function(x, y, type, ...) {
panel.xYplot(x, y, type=p,...)
lpoints(x, y, pch=16, col=white, cex=2)
panel.superpose(x, y, type=p, ...)
},
par.settings = sp,
auto.key=list(columns=1, x=0.7, y=0.8, corner = c(0,0))
)
Sander Oom wrote:
An off list response from Mat Soukop (thanks Mat!!) provides an even
more elegant solution (see code below)! I have included the original
code, so people can decide whether to plot in a single panel or in
multiple panels. Now we have a fully functional workaround to get
plotmeans{gplots} for multiple factors using lattice! Great!
library(Hmisc)
library(lattice)
tmp -
structure(list(Position = structure(as.integer(c(1, 2, 1, 2,
1, 2, 1, 2)), .Label = c(Inside, Outside), class = factor),
AltGeo = structure(as.integer(c(1, 1, 2, 2, 3, 3, 4, 4)), .Label =
c(Basalt-High,
Basalt-Low, Quartz-High, Quartz-Low), class = factor),
Sodium = c(27.3, 26.9, 25, 18.1,
4.67, 5.56, 10.7, 5.67
), SD = c(5.3851648071345, 2.42097317438899, 20.1618451536560,
15.2679766541317, 5.45435605731786, 8.09492296305393, 10.6183802907976,
8.06225774829855), Nobs = c(9, 9, 9, 9, 9, 9, 9, 9), Lower =
c(25.5382783976218,
26.0818978307592, 18.2793849487813, 13.0217855597339, 2.84854798089405,
2.85724790120425, 7.12720656973412, 2.97924741723382), Upper =
c(29.1283882690448,
27.6958799470186, 31.7206150512187, 23.2004366624884, 6.48478535243929,
8.25386320990686, 14.2061267635992, 8.35408591609952)), .Names =
c(Position,
AltGeo, Sodium, SD, Nobs, Lower, Upper), row.names = c(1,
2, 3, 4, 5, 6, 7, 8), class = data.frame)
tmp$PosNum - unclass(tmp$Position)
tmp
(labs - unique(tmp$Position))
# plot factor levels in seperate panels
xYplot(Cbind(Sodium,Lower,Upper) ~ PosNum|AltGeo, data=tmp, nx=FALSE,
xlim=c(0.5,2.5),
ylim=c(min(tmp$Lower)-1,max(tmp$Upper)+1),
scales = list(x = list(at=seq(1, 2, by=1), labels = labs)),
xlab=Position, ylab=Sodium
)
new.back - trellis.par.get(background)
new.back$col - white
newcol - trellis.par.get(superpose.symbol)
newcol$col - c('green4','blue','red','black')
newcol$pch - c(16,1,4,8)
new.line - trellis.par.get(box.rectangle)
new.line$col - 'black'
trellis.par.set(background, new.back)
trellis.par.set(superpose.symbol, newcol)
trellis.par.set(box.rectangle, new.line)
# Plot factor levels in one graph
tmp$xvar - rep(1:4, each=2)+rep(c(-.05,.05), 4)
xYplot(Cbind(Sodium,Lower,Upper) ~ xvar, groups=Position, data=tmp,
scales=list(y='free',x=list(at=1:4,
labels=levels(tmp$AltGeo))),
xlab='AltGeo', xlim=c(.5, 4.5),
key=list(points=Rows(trellis.par.get(superpose.symbol),1:2),
text=list(lab =as.character(levels(tmp$Position)),
col=trellis.par.get(superpose.symbol)$col[1:2]),
columns=2, cex=1, title=Position,
cex.title=1.1))
Sander Oom wrote:
Problem solved!
I was so focused on reproducing
Re: [R] Problem with data frame when using xYplot?
Hi Deepayan!
Deepayan Sarkar wrote:
On Friday 13 May 2005 08:07 am, Sander Oom wrote:
An off list response from Mat Soukop (thanks Mat!!) provides an even
more elegant solution (see code below)! I have included the original
code, so people can decide whether to plot in a single panel or in
multiple panels. Now we have a fully functional workaround to get
plotmeans{gplots} for multiple factors using lattice! Great!
Just out of curiousity, does replacing the 'xlim' and 'scales' arguments above
by
xlim = levels(tmp$Position)
do the same thing? It should with xyplot (which also allows the x variable to
be a factor), but xYplot may be bypassing that.
You mean: xlim = levels(tmp$AltGeo)yes it does!? No clue how one
would ever get comfortable with all these options!
library(Hmisc)
library(lattice)
ltheme - canonical.theme(color = FALSE) ## in-built BW theme
ltheme$strip.background$col - transparent ## change strip bg
lattice.options(default.theme = ltheme) ## set as default
tmp -
structure(list(Position = structure(as.integer(c(1, 2, 1, 2,
1, 2, 1, 2)), .Label = c(Inside, Outside), class = factor),
AltGeo = structure(as.integer(c(1, 1, 2, 2, 3, 3, 4, 4)), .Label =
c(Basalt-High,
Basalt-Low, Quartz-High, Quartz-Low), class = factor),
Sodium = c(27.3, 26.9, 25, 18.1,
4.67, 5.56, 10.7, 5.67
), SD = c(5.3851648071345, 2.42097317438899, 20.1618451536560,
15.2679766541317, 5.45435605731786, 8.09492296305393, 10.6183802907976,
8.06225774829855), Nobs = c(9, 9, 9, 9, 9, 9, 9, 9), Lower =
c(25.5382783976218,
26.0818978307592, 18.2793849487813, 13.0217855597339, 2.84854798089405,
2.85724790120425, 7.12720656973412, 2.97924741723382), Upper =
c(29.1283882690448,
27.6958799470186, 31.7206150512187, 23.2004366624884, 6.48478535243929,
8.25386320990686, 14.2061267635992, 8.35408591609952)), .Names =
c(Position,
AltGeo, Sodium, SD, Nobs, Lower, Upper), row.names = c(1,
2, 3, 4, 5, 6, 7, 8), class = data.frame)
tmp$xvar - rep(1:4, each=2)+rep(c(-0.1,0.1), 4)
tmp
sp - list(superpose.symbol = list(pch = c(16,1), cex = 1))
xYplot(Cbind(Sodium,Lower,Upper) ~ xvar, groups=Position, data=tmp,
xlim = levels(tmp$AltGeo),
ylim=c(min(tmp$Lower)-1,max(tmp$Upper)+1),
xlab='AltGeo', ylab='Sodium',
panel = function(x, y, type, ...) {
panel.xYplot(x, y, type=p,...)
lpoints(x, y, pch=16, col=white, cex=2)
panel.superpose(x, y, type=p, ...)
},
par.settings = sp,
auto.key=list(columns=1, x=0.7, y=0.8, corner = c(0,0))
)
Deepayan
__
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
--
Dr. Sander P. Oom
Animal, Plant and Environmental Sciences,
University of the Witwatersrand
Private Bag 3, Wits 2050, South Africa
Tel (work) +27 (0)11 717 64 04
Tel (home) +27 (0)18 297 44 51
Fax +27 (0)18 299 24 64
Email [EMAIL PROTECTED]
Web www.oomvanlieshout.net/sander
__
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
Re: [R] Problem with data frame when using xYplot?
Sorry Deepayan, I forgot that the code moved on while you send your
reply. Below the simplified version using your suggestion and this time
based on the generic data used in the xYplot manual!
Maybe this example can be included in the manual, so next time people
will find the answer there. Or better still, you could provide a high
level function in 'lattice'. ;-)
dfr - expand.grid(month=1:12, continent=c('Europe','USA'),
sex=c('female','male'))
set.seed(1)
dfr$y - dfr$month/10 + 1*(dfr$sex=='female') +
2*(dfr$continent=='Europe') + runif(48,-.15,.15)
dfr
dfs - summarize(dfr$y, llist(dfr$continent,dfr$sex),
smean.cl.normal)
labs - unique(dfs$continent)
colnames(dfs) - list(continent,sex,y,Lower,Upper)
dfs$sexnum - unclass(dfs$sex)
dfs
xYplot(Cbind(y,Lower,Upper) ~ sexnum|continent, data=dfs, nx=FALSE,
xlim=levels(dfs$sex),
ylim=c(min(dfs$Lower)-1,max(dfs$Upper)+1),
)
dfs$xvar - rep(1:2, each=2)+rep(c(-0.1,0.1), 2)
dfs
sp - list(superpose.symbol = list(pch = c(16,1), cex = 1),
superpose.line = list(col = grey, lty = 1))
xYplot(Cbind(y,Lower,Upper) ~ xvar, groups=sex, data=dfs,
xlim= levels(dfs$continent),
ylim= c(min(dfs$Lower)-1,max(dfs$Upper)+1),
xlab=Continent,
panel=function(x, y, type, ...) {
panel.xYplot(x, y, type=p,...)
lpoints(x, y, pch=16, col=white, cex=2)
panel.superpose(x, y, type=p, ...)
},
par.settings= sp,
auto.key= list(columns=1, x=0.7, y=0.8, corner = c(0,0))
)
Deepayan Sarkar wrote:
On Friday 13 May 2005 10:36 am, Sander Oom wrote:
Hi Deepayan!
Deepayan Sarkar wrote:
On Friday 13 May 2005 08:07 am, Sander Oom wrote:
An off list response from Mat Soukop (thanks Mat!!) provides an even
more elegant solution (see code below)! I have included the original
code, so people can decide whether to plot in a single panel or in
multiple panels. Now we have a fully functional workaround to get
plotmeans{gplots} for multiple factors using lattice! Great!
Just out of curiousity, does replacing the 'xlim' and 'scales' arguments
above by
xlim = levels(tmp$Position)
do the same thing? It should with xyplot (which also allows the x
variable to be a factor), but xYplot may be bypassing that.
You mean: xlim = levels(tmp$AltGeo)yes it does!?
No, I meant exactly what I wrote, and my comment followed this piece of code
(which you have deleted from your reply):
-
tmp$PosNum - unclass(tmp$Position)
tmp
(labs - unique(tmp$Position))
# plot factor levels in seperate panels
xYplot(Cbind(Sodium,Lower,Upper) ~ PosNum|AltGeo, data=tmp, nx=FALSE,
xlim=c(0.5,2.5),
ylim=c(min(tmp$Lower)-1,max(tmp$Upper)+1),
scales = list(x = list(at=seq(1, 2, by=1), labels = labs)),
xlab=Position, ylab=Sodium
)
--
No clue how one would ever get comfortable with all these options!
By reading the manual, of course :-)
Deepayan
__
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
--
Dr. Sander P. Oom
Animal, Plant and Environmental Sciences,
University of the Witwatersrand
Private Bag 3, Wits 2050, South Africa
Tel (work) +27 (0)11 717 64 04
Tel (home) +27 (0)18 297 44 51
Fax +27 (0)18 299 24 64
Email [EMAIL PROTECTED]
Web www.oomvanlieshout.net/sander
__
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
Re: [R] Conflict between xtable and Hmisc when using Sweave?
Dear Frank,
I have a Sweave document in which I export anova (aov) tables to Latex
and calculate some summary statistics with summarize{Hmisc} for a graph
(as in the example below).
I currently use the following code for the aov tables:
results=tex=
tmp - datGrassHC[datGrassHC$Loc 0 datGrassHC$Loc 9 ,]
tmpAov - aov(Height~Geology*Altitude*Origin*BinInOut , data=tmp)
tmpTable - xtable (tmpAov ,
caption=ANOVA table for vegetation height.,
label=tab:AnovaHeight
)
print.xtable(tmpTable, type=latex, floating=TRUE,
table.placement=ht, caption.placement=top,
latex.environments=c(center))
)
@
I used xtables, because it has a working aov example. I would be happy
to use an alternative if I knew how! Would you have sample code to
illustrate how to export an aov table to Latex using latex{Hmisc}.
Thanks very much for your help,
Sander.
Frank E Harrell Jr wrote:
Sander Oom wrote:
Dear R users,
The Sweave code below runs fine, as it is. However, an error occurs
when the line 'library(xtable)' is uncommented:
Error: chunk 1
Error in label-(`*tmp*`, value = month) :
no applicable method for label-
Is anybody aware of this and knows a workaround?
Thanks,
Sander.
***
\documentclass[a4paper]{article}
\title{Sweave Test for summarize}
\author{Sander Oom}
\usepackage{a4wide}
\begin{document}
\maketitle
\begin{figure}[ht]
\begin{center}
fig=TRUE,echo=FALSE=
# library(xtable)
library(Hmisc)
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)
print(xYplot(Cbind(y,Lower,Upper) ~ month, groups=year, data=s,
keys='lines', method='alt', type='b'))
@
\end{center}
\end{figure}
\end{document}
version
_
platform i686-pc-linux-gnu
arch i686
os linux-gnu
system i686, linux-gnu
status
major2
minor1.0
year 2005
month04
day 18
language R
I feel this is an xtable problem because Hmisc has being using label and
label- since 1991.
Frank
--
Dr. Sander P. Oom
Animal, Plant and Environmental Sciences,
University of the Witwatersrand
Private Bag 3, Wits 2050, South Africa
Tel (work) +27 (0)11 717 64 04
Tel (home) +27 (0)18 297 44 51
Fax +27 (0)18 299 24 64
Email [EMAIL PROTECTED]
Web www.oomvanlieshout.net/sander
__
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
[R] Problem with data frame when using xYplot?
Dear all, I am trying to plot means and error bars using xYplot, but I get an error message from xYplot which I can not figure out: Error in Summary.factor(..., na.rm = na.rm) : range not meaningful for factors The data frame (tmpNa) was created using aggregate. I have used dump to created the code below, which generates the same error. Can anybody tell me what is wrong with the data frame? Thanks in advance, Sander. library(Hmisc) tmpNa - structure(list(Position = structure(as.integer(c(1, 2, 1, 2, 1, 2, 1, 2)), .Label = c(Inside, Outside), class = factor), AltGeo = structure(as.integer(c(1, 1, 2, 2, 3, 3, 4, 4)), .Label = c(Basalt-High, Basalt-Low, Quartz-High, Quartz-Low), class = factor), Sodium = c(27.3, 26.9, 25, 18.1, 4.67, 5.56, 10.7, 5.67 ), SD = c(5.3851648071345, 2.42097317438899, 20.1618451536560, 15.2679766541317, 5.45435605731786, 8.09492296305393, 10.6183802907976, 8.06225774829855), Nobs = c(9, 9, 9, 9, 9, 9, 9, 9), Lower = c(25.5382783976218, 26.0818978307592, 18.2793849487813, 13.0217855597339, 2.84854798089405, 2.85724790120425, 7.12720656973412, 2.97924741723382), Upper = c(29.1283882690448, 27.6958799470186, 31.7206150512187, 23.2004366624884, 6.48478535243929, 8.25386320990686, 14.2061267635992, 8.35408591609952)), .Names = c(Position, AltGeo, Sodium, SD, Nobs, Lower, Upper), row.names = c(1, 2, 3, 4, 5, 6, 7, 8), class = data.frame) xYplot(Cbind(Sodium,Lower,Upper) ~ AltGeo, groups=Position, data=tmpNa) version _ platform i686-pc-linux-gnu arch i686 os linux-gnu system i686, linux-gnu status major2 minor1.0 year 2005 month04 day 18 language R -- Dr. Sander P. Oom Animal, Plant and Environmental Sciences, University of the Witwatersrand Private Bag 3, Wits 2050, South Africa Tel (work) +27 (0)11 717 64 04 Tel (home) +27 (0)18 297 44 51 Fax +27 (0)18 299 24 64 Email [EMAIL PROTECTED] Web www.oomvanlieshout.net/sander __ 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
Re: [R] Plotting means and confidence intervals by group factor using lattice graphics?
Thanks for your help Deepayan!
I noticed the other thread before, but I was really looking for a dot
plot with error bars. This looks better indeed, as you state yourself in
that thread.
The package Hmisc gives the solution through Cbind and xYplot! I will
prepare an example graph for the gallery.
Thanks again,
Sander.
library(Hmisc)
# Examples of plotting raw data
dfr - expand.grid(month=1:12, continent=c('Europe','USA'),
sex=c('female','male'))
set.seed(1)
dfr - upData(dfr,
y=month/10 + 1*(sex=='female') + 2*(continent=='Europe') +
runif(48,-.15,.15),
lower=y - runif(48,.05,.15),
upper=y + runif(48,.05,.15))
xYplot(Cbind(y,lower,upper) ~ month,subset=sex=='male' continent=='USA',
data=dfr)
xYplot(Cbind(y,lower,upper) ~ month|continent, subset=sex=='male',data=dfr)
xYplot(Cbind(y,lower,upper) ~ month|continent, groups=sex, data=dfr); Key()
# add ,label.curves=FALSE to suppress use of labcurve to label curves where
# farthest apart
Deepayan Sarkar wrote:
On Wednesday 04 May 2005 10:30, Sander Oom wrote:
Dear R graphics gurus,
Another question about lattice graphics. This time I would like to plot
means and confidence intervals by group factor in a lattice graph. I can
not find any working lattice examples. Maybe a custom panel function is
the answer, but that is a bit beyond me for now.
There's an example in this thread:
http://finzi.psych.upenn.edu/R/Rhelp02a/archive/50299.html
There might be useful tools in the Hmisc package as well.
The individual plots within the lattice graph could look like this:
# Example with confidence intervals and grid
hh - t(VADeaths)[, 5:1]
mybarcol - gray20
ci.l - hh * 0.85
ci.u - hh * 1.15
mp - barplot2(hh, beside = TRUE,
col = c(lightblue, mistyrose,
lightcyan, lavender),
legend = colnames(VADeaths), ylim = c(0, 100),
main = Death Rates in Virginia, font.main = 4,
sub = Faked 95 percent error bars, col.sub = mybarcol,
cex.names = 1.5, plot.ci = TRUE, ci.l = ci.l, ci.u = ci.u,
plot.grid = TRUE)
This gives me
Error: couldn't find function barplot2
Maybe you missed a library() call?
Deepayan
mtext(side = 1, at = colMeans(mp), line = -2,
text = paste(Mean, formatC(colMeans(hh))), col = red)
box()
Or like this:
data(state)
plotmeans(state.area ~ state.region)
Both plotmeans and barplot2 give interesting options such as printing of
nobs, among other things. In case of a barplot, there should be an
option to plot the confidence intervals in one direction only (up) as to
avoid interference with any black and white shading. The plotMeans
function provides a useful option error.bars (se, sd, conf.int,
none).
The following test data is still useful:
tmp - expand.grid(geology = c(Sand,Clay,Silt,Rock),
species =
c(ArisDiff,BracSera,CynDact,ElioMuti,EragCurS,EragPseu),
dist = seq(1,9,1) )
tmp$height - rnorm(216)
For instance plotting height versus dist by geology.
Any help very welcome!
Cheers,
Sander.
PS Of course the resulting graph will go to the R graph gallery!
__
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
--
Dr. Sander P. Oom
Animal, Plant and Environmental Sciences,
University of the Witwatersrand
Private Bag 3, Wits 2050, South Africa
Tel (work) +27 (0)11 717 64 04
Tel (home) +27 (0)18 297 44 51
Fax +27 (0)18 299 24 64
Email [EMAIL PROTECTED]
Web www.oomvanlieshout.net/sander
__
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
[R] Plotting means and confidence intervals by group factor using lattice graphics?
Dear R graphics gurus, Another question about lattice graphics. This time I would like to plot means and confidence intervals by group factor in a lattice graph. I can not find any working lattice examples. Maybe a custom panel function is the answer, but that is a bit beyond me for now. The individual plots within the lattice graph could look like this: # Example with confidence intervals and grid hh - t(VADeaths)[, 5:1] mybarcol - gray20 ci.l - hh * 0.85 ci.u - hh * 1.15 mp - barplot2(hh, beside = TRUE, col = c(lightblue, mistyrose, lightcyan, lavender), legend = colnames(VADeaths), ylim = c(0, 100), main = Death Rates in Virginia, font.main = 4, sub = Faked 95 percent error bars, col.sub = mybarcol, cex.names = 1.5, plot.ci = TRUE, ci.l = ci.l, ci.u = ci.u, plot.grid = TRUE) mtext(side = 1, at = colMeans(mp), line = -2, text = paste(Mean, formatC(colMeans(hh))), col = red) box() Or like this: data(state) plotmeans(state.area ~ state.region) Both plotmeans and barplot2 give interesting options such as printing of nobs, among other things. In case of a barplot, there should be an option to plot the confidence intervals in one direction only (up) as to avoid interference with any black and white shading. The plotMeans function provides a useful option error.bars (se, sd, conf.int, none). The following test data is still useful: tmp - expand.grid(geology = c(Sand,Clay,Silt,Rock), species = c(ArisDiff,BracSera,CynDact,ElioMuti,EragCurS,EragPseu), dist = seq(1,9,1) ) tmp$height - rnorm(216) For instance plotting height versus dist by geology. Any help very welcome! Cheers, Sander. PS Of course the resulting graph will go to the R graph gallery! -- Dr. Sander P. Oom Animal, Plant and Environmental Sciences, University of the Witwatersrand Private Bag 3, Wits 2050, South Africa Tel (work) +27 (0)11 717 64 04 Tel (home) +27 (0)18 297 44 51 Fax +27 (0)18 299 24 64 Email [EMAIL PROTECTED] Web www.oomvanlieshout.net/sander __ 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
[R] Help needed with lattice graph!
Dear R users, If I manage to sort out this graph, it is certainly a candidate for the new R graph gallery (http://addictedtor.free.fr/graphiques/displayGallery.php)! I created the following lattice graph: library(lattice) tmp - expand.grid(geology = c(Sand,Clay,Silt,Rock), species = c(ArisDiff,BracSera,CynDact,ElioMuti,EragCurS,EragPseu), dist = seq(1,9,1) ) tmp$height - rnorm(216) sps - trellis.par.get(superpose.symbol) sps$pch - 1:6 trellis.par.set(superpose.symbol, sps) xyplot( height ~ dist | geology, data = tmp, groups = species, type = b, cex = 1.2, layout = c(2,2), lines = list(col=grey), key = list(columns = 2, type = b, cex = 1.2, text = list(paste(unique(tmp$species))), points = Rows(sps, 1:6) ) ) However, for once, the R defaults are not to my liking. I plot the graph to postscript and the result is less then optimal. I would like to plot the point symbols in black and white, both in the graphs and the key. I would like the lines to be a single style (grey or a light dash) and preferably the lines do not go through the symbols (like figure 4.11 in the MASS book). I have tried many, many options, but results varied from wrong symbols to wrong things plotted. Splitting the lines and points over different panels seems the way, but I can not make it work. Your help is much appreciated! This graph and the resulting black and white graph will be posted on the R graph gallery. Thanks, Sander. -- Dr. Sander P. Oom Animal, Plant and Environmental Sciences, University of the Witwatersrand Private Bag 3, Wits 2050, South Africa Tel (work) +27 (0)11 717 64 04 Tel (home) +27 (0)18 297 44 51 Fax +27 (0)18 299 24 64 Email [EMAIL PROTECTED] Web www.oomvanlieshout.net/sander __ 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
[R] Landscape indeces analysis methods as an R package!?
Dear Barry, Thanks for your stimulating reply. As you see I have send a CC of this reply to the R mailing list. The R mailing list will be the best place to send your queries. You can subscribe here: http://www.r-project.org/mail.html I would think there is an interest in the provision of landscape indices analysis methods in R. There are several packages for landscape indices analysis available (Fragstats, APACK, IAN), but all use their own user and data interfaces. Only the r.le tools are available in a more general software: GRASS GIS. It would be usefull to have the functionality within the powerful data analysis tool R. Of course providing an R package will remove most of the user and data interface requirements as they will be taken care off by R. Only the core landscape indices methods need to be included. You will have to consult the R mailing list on how to translate your code into an R package! The spatial aspects of R are further developed at the moment and spatial packages become more and more powerful. Look here for most recent development: http://www.r-project.org/Rgeo R can currently read (relevant R package name in brackets: ASCII grid (with custom function or Adehabitat), shape files (maptools and shapefiles), GRASS files (grass) and ArcInfo files (RArcInfo), many others via external gdal (rgdal). I would say plenty of opportunities! Hope this gets you under way! Cheers, Sander. DeZonia, Barry wrote: Hi Sander, APACK is no longer being developed but is still available for download. There isn't a Linux version. APACK has been superceded by IAN (http://landscape.forest.wisc.edu/projects/ian/). It is written in Ruby and runs on Windows platforms. A Linux port could be done with a little effort. You are only the second person to talk about making an R package. I've considered this and may yet do so. I've been studying what to do and how for a little while now. If you wanted to provide input as to what R users' needs are I would be appreciative. One question I have is how do users currently get raster data loaded into R (as a matrix I suppose but what file formats are supported)? I've seen the pixmap library but the formats are not broadly supported. -Original Message- From: Sander Oom [mailto:[EMAIL PROTECTED] Sent: Wednesday, March 16, 2005 3:49 AM To: [EMAIL PROTECTED] Subject: Apack information? Dear Apack developers, Could you inform me about the current status of APACK? Is a version available for Linux? Would it be possible to translate APACK into an R package? Looking forward to more information. Yours, Sander Oom. __ 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
Re: [R] from long/lat to UTM
Hi Yyan, The proj4R package by Roger Bivand will allow you to project data in many ways and directions. http://spatial.nhh.no/R/Devel/proj4R-pkg.pdf It uses the proj libraries from: http://www.remotesensing.org/proj/ Not sure where you would derive the time zone! Good luck, Sander. yyan liu wrote: Hi: Is there any function in R which can convert the long/lat to UTM(Universal Transverse Mercator)? There are quite a few converters on Internet. However, the interface is designed as input-output which I can not convert lots of locations at the same time. Another question is whether there is a function in R which can tell the time zone from the location's lat/long? Thank you! liu __ 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 -- Dr. Sander P. Oom Animal, Plant and Environmental Sciences, University of the Witwatersrand Private Bag 3, Wits 2050, South Africa Tel (work) +27 (0)11 717 64 04 Tel (home) +27 (0)18 297 44 51 Fax +27 (0)18 299 24 64 Email [EMAIL PROTECTED] Web www.oomvanlieshout.net/sander __ 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
Re: [R] zeichnung eines wahrscheinlichkeitsnetztes
Hi Ernst, Pretty soon some one is going to: 1) tell you that the mailing list is in English, 2) that you should provide more details of your data and 2) that your message should show signs that you have read the posting guide and have searched the help documentation in a first attempt to find the answer to your question. Maybe you can reply to you own question, considering the above. Then people will surely be very helpful! Good luck, Sander. ernst hoffer wrote: ich hätte ein problem: im rahmen meiner projektarbeit stehe ich vor der aufgabe, einen beliebigen datensatzt einzulesen, die empirische verteilungsfunktion des datensatztes zu zeichnen, aber das ganze so, dass die y-achse nach der quantilsfunktion der exponetialverteilung mit parameter=1 transformiert werden soll! leider steh ich da jetzt wirklich an! würde mich serh freuen, wenn mir jemand weiterhelfen könnte, mit herzlichem dank, reka [[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 -- - Dr. Sander P. Oom Animal, Plant and Environmental Sciences University of the Witwatersrand Private Bag 3 Wits 2050 South Africa Tel (work) +27 (0)11 717 64 04 Tel (home) +27 (0)18 297 44 51 Fax +27 (0)18 299 24 64 Email [EMAIL PROTECTED] Web www.oomvanlieshout.net/sander __ 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
Re: [R] Hosting a R Graph Gallery?
Thanks for the responses and offer! Romain, good luck with your exams first of all! Graphics in R base and R contributed is a good start indeed. I thought of a wiki as well as it will require less maintenance from the host. Custom html could be powerful, but will require more input from host! With build-in search engines, indexing etc., a wiki environment could provide all functionality for an easy to navigate gallery! To stimulate development and support of SVG functionality, we might want to offer people the option of submitting graphics in both png and svg! The work by Jake looks very promising: http://www.darkridge.com/~jake/RSvg/ The applet works well with the examples! Cheers, Sander. Robert Cunningham wrote: I too have often though a R-gallery would be useful. It seems to me that a Wiki-style page with a database backend would be the best bet. It also seems to be that the best place to start is a complete image gallery produced from all the examples in R base, then in packages in CRAN. In this context the graphicsQC package (http://www.stat.auckland.ac.nz/~paul/R/graphicsQC_0.4.tar.g) of Paul Murrell seems useful. Cheers, Robert Cunningham Romain Francois [EMAIL PROTECTED] writes: Hello Sander, That's a good idea and i am up to it. Right now i am in an exam period, so it's not really the better time, give me a couple of weeks and i will come up with a specific format of R files to submit to me that i could post-process to generate html documents. To my mind, those html files should show : - the plot itself + Submitter(s) - web page - email (eventually protected, I don't know how to do it) - Bibliographic references - Required R packages + Commentaries - in english - and in any other languages I'm open to any suggestion. Romain. Le 18.02.2005 14:33, Sander Oom a écrit : Dear R users, Following some of the recent questions and discussions about the R plotting abilities, it occurred to me again that it would be very valuable to have an R graph gallery. Eric Lecoutre made a very nice example in: http://www.stat.ucl.ac.be/ISpersonnel/lecoutre/stats/fichiers/_gallery.pdf It would be very useful to many beginners, but probably also advanced users of R, to have an overview of R graph types with graphical examples and associated R code. In order to facilitate the evolution of a large gallery, some sort of wiki environment might be most suitable, thus providing access to all users, but with limited maintenance costs for the provider. Do others agree this could be a valuable resource? Would anybody have the resources to host such an R graph gallery? Yours, Sander Oom. -- Romain FRANCOIS : [EMAIL PROTECTED] page web : http://addictedtor.free.fr/ (en construction) 06 18 39 14 69 / 01 46 80 65 60 ___ Etudiant en 3eme année Institut de Statistique de l'Université de Paris (ISUP) Filière Industrie et Services http://www.isup.cicrp.jussieu.fr/ __ 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 __ 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 -- - Dr. Sander P. Oom Animal, Plant and Environmental Sciences University of the Witwatersrand Private Bag 3 Wits 2050 South Africa Tel (work) +27 (0)11 717 64 04 Tel (home) +27 (0)18 297 44 51 Fax +27 (0)18 299 24 64 Email [EMAIL PROTECTED] Web www.oomvanlieshout.net/sander __ 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
[R] Repeating grey scale in graph?
Dear R users, Could somebody tell me why the grey color ramp is repeated in this graph, eventhough the ramp values go from 0 to 1? I must be missing something obvious, but I can not see it! z - c(0.064329041,0.117243316,0.161565116,0.19923015,0.231642175,0.259835539,0.284571226, 0.038507288,0.094184749,0.140959431,0.180803984,0.215159105,0.245096084,0.271412845, 0.00775022,0.066198255,0.115433207,0.157494219,0.193836765,0.225569076,0.253518629, -0.02820814,0.032958752,0.084661362,0.128946221,0.167320522,0.200892494,0.230504392, -0.07003273,-0.005814512,0.048304039,0.094805358,0.135196637,0.170630435,0.201956395, -0.117878701,-0.050461393,0.005991829,0.054672666,0.097103088,0.134398711,0.167423957) x - c(0,1,2,3,4,5) y - c(50, 100, 150, 200, 250, 300, 350) z - matrix(z, nrow=length(x), ncol=length(y), byrow=TRUE) #persp(x, y, z, theta = 30, phi = 30, expand = 0.5, # box= TRUE, axes= TRUE, ticktype = detailed, main=Title of plot) hgt - (z - min(z))/ (max(z) - min(z)) z hgt cols - grey(hgt) persp(x, y, z, col = cols, theta = 30, phi = 30, expand = 0.5, box= TRUE, axes= TRUE, ticktype = detailed, main=Title of plot) Thanks, Sander. version _ platform i386-pc-mingw32 arch i386 os mingw32 system i386, mingw32 status major2 minor0.1 year 2004 month11 day 15 language R -- - Dr. Sander P. Oom Animal, Plant and Environmental Sciences University of the Witwatersrand Private Bag 3 Wits 2050 South Africa Tel (work) +27 (0)11 717 64 04 Tel (home) +27 (0)18 297 44 51 Fax +27 (0)18 299 24 64 Email [EMAIL PROTECTED] Web www.oomvanlieshout.net/sander __ 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
Re: [R] Repeating grey scale in graph?
Thanks Peter! Of course I only have (nx-1)(ny-1) facets in a x*y plot! The help page line: ... col the color(s) of the surface facets. Transparent colours are ignored. This is recycled to the (nx-1)(ny-1) facets. ... just did not ring a bell. In fact, it is still not clear to me why it recycles the ramp even though it has a surplus of colours (grey levels)! Why not just ignore the surplus colours? Anyway it works, Sander. Peter Dalgaard wrote: Sander Oom [EMAIL PROTECTED] writes: Dear R users, Could somebody tell me why the grey color ramp is repeated in this graph, eventhough the ramp values go from 0 to 1? I must be missing something obvious, but I can not see it! z - c(0.064329041,0.117243316,0.161565116,0.19923015,0.231642175,0.259835539,0.284571226, 0.038507288,0.094184749,0.140959431,0.180803984,0.215159105,0.245096084,0.271412845, 0.00775022,0.066198255,0.115433207,0.157494219,0.193836765,0.225569076,0.253518629, -0.02820814,0.032958752,0.084661362,0.128946221,0.167320522,0.200892494,0.230504392, -0.07003273,-0.005814512,0.048304039,0.094805358,0.135196637,0.170630435,0.201956395, -0.117878701,-0.050461393,0.005991829,0.054672666,0.097103088,0.134398711,0.167423957) x - c(0,1,2,3,4,5) y - c(50, 100, 150, 200, 250, 300, 350) z - matrix(z, nrow=length(x), ncol=length(y), byrow=TRUE) #persp(x, y, z, theta = 30, phi = 30, expand = 0.5, # box= TRUE, axes= TRUE, ticktype = detailed, main=Title of plot) hgt - (z - min(z))/ (max(z) - min(z)) z hgt cols - grey(hgt) persp(x, y, z, col = cols, theta = 30, phi = 30, expand = 0.5, box= TRUE, axes= TRUE, ticktype = detailed, main=Title of plot) You have 30 facets and 42 colour values. Try it with cols - grey(hgt[-1,-1]) -- - Dr. Sander P. Oom Animal, Plant and Environmental Sciences University of the Witwatersrand Private Bag 3 Wits 2050 South Africa Tel (work) +27 (0)11 717 64 04 Tel (home) +27 (0)18 297 44 51 Fax +27 (0)18 299 24 64 Email [EMAIL PROTECTED] Web www.oomvanlieshout.net/sander __ 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
Re: [R] Repeating grey scale in graph?
Aaaah...the inner workings of R! Now I also see why the colours are not only repeated, but also 'wrongly' allocated to the facets! Very clear example! Indeed a warning or error would have been more helpful! Cheers, Sander. PS: I hope that after all this, I can still convince the creator of the original data that it is a good idea to plot his graphs in R instead of excel. ;-) Duncan Murdoch wrote: On Wed, 16 Feb 2005 15:44:00 +0200, Sander Oom [EMAIL PROTECTED] wrote : Thanks Peter! Of course I only have (nx-1)(ny-1) facets in a x*y plot! The help page line: ... col the color(s) of the surface facets. Transparent colours are ignored. This is recycled to the (nx-1)(ny-1) facets. ... just did not ring a bell. In fact, it is still not clear to me why it recycles the ramp even though it has a surplus of colours (grey levels)! Why not just ignore the surplus colours? Your z array is 6 by 7. Your cols will be mapped to a 5 by 6 array. They don't look like an array, because the grey() function stripped off the dimension attribute. The problem is that if you pass the entries from a 6 by 7 array to something that expects the entries from a 5 by 6 array, you get things in the wrong order. You see the same effect here: rownum - as.vector(row(matrix(NA, 6, 7))) matrix(rownum, 6, 7) [,1] [,2] [,3] [,4] [,5] [,6] [,7] [1,]1111111 [2,]2222222 [3,]3333333 [4,]4444444 [5,]5555555 [6,]6666666 matrix(rownum, 5, 6) [,1] [,2] [,3] [,4] [,5] [,6] [1,]165432 [2,]216543 [3,]321654 [4,]432165 [5,]543216 Warning message: data length [42] is not a sub-multiple or multiple of the number of rows [5] in matrix except that in this case you get a warning about the wrong length; persp doesn't give you the warning. Maybe it should? Duncan Murdoch -- - Dr. Sander P. Oom Animal, Plant and Environmental Sciences University of the Witwatersrand Private Bag 3 Wits 2050 South Africa Tel (work) +27 (0)11 717 64 04 Tel (home) +27 (0)18 297 44 51 Fax +27 (0)18 299 24 64 Email [EMAIL PROTECTED] Web www.oomvanlieshout.net/sander __ 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
Re: [R] Time line plot in R?
Jim,
Brilliant! Thought someone might have figured it out already. Now we
just need a gallery to show off this graph!
One little thing: the par(mar=6,6,4,2) gives an error:
'Error in par(args) : parameter mar has the wrong length'
Any suggestions?
Code below includes fake labels for testing.
Cheers,
Sander.
# some fake data please, maestro
fakedata-sample(1:10,10)
# leave a bit of extra space beneath and to the left of the plot
par(mar=6,6,4,2)
# this function will probably end up in the plotrix package
time.line-function(x,y,at=NULL,labels=TRUE,tlticks=NULL,...) {
if(missing(x) missing (y))
stop(Usage: time.line(x,y,at=NULL,labels=TRUE))
plotrange-par(usr)
# get the graphical parameters
oldpar-par(no.readonly=TRUE)
# turn off clipping
par(xpd=TRUE)
if(missing(x)) {
# it's a horizontal line
segments(plotrange[1],y,plotrange[2],y,...)
ticklen-(plotrange[4]-plotrange[3])*0.02
if(!is.null(tlticks))
segments(tlticks,y+ticklen,tlticks,y-ticklen,...)
mwidth-strwidth(M)
# blank out the line where labels will appear
rect(at-mwidth,y-ticklen,at+mwidth,y+ticklen,col=white,border=FALSE)
# rotate the text
par(srt=270)
# draw the labels
text(at,y,labels,...)
}
if(missing(y)) {
# it's a vertical line
# draw the line
segments(x,plotrange[3],x,plotrange[4],...)
ticklen-(plotrange[2]-plotrange[1])*0.02
if(!is.null(tlticks))
segments(x+ticklen,tlticks,x-ticklen,tlticks,...)
mheight-strheight(M)
# blank out the line where labels will appear
rect(x-ticklen,at-mheight,x+ticklen,at+mheight,col=white,border=FALSE)
# draw the labels
text(x,at,labels,...)
}
# restore the parameters
par(oldpar)
}
# some fake labels
eventT-c(2.5, 4, 7, 8.5)
eventD-c(first event,second event,third event,fourth event)
tl.labels-data.frame(eventT,eventD)
plot(fakedata,xlab=)
# display a horizontal time line
time.line(x=-1,at=tl.labels[[1]],labels=tl.labels[[2]],col=black)
# now a vertical one
time.line(y=-1,at=tl.labels[[1]],labels=tl.labels[[2]],col=black)
version
_
platform i386-pc-mingw32
arch i386
os mingw32
system i386, mingw32
status
major2
minor0.1
year 2004
month11
day 15
language R
Jim Lemon wrote:
Sander Oom wrote:
Dear R users,
In order to illustrate the possible effects of events on variables
plotted against time, I would like plot a time line of events along side
the plot of the variables.
The x-axis should be some time unit; the y-axis should be the variable
of interest; the time line should be plotted below the graph along the
same x-axis scale.
As I have many variables and different events lists, I would like to
write a script to read the events from a file and plot them together
with the other plot.
The time line should look something like this:
http://www.oslis.k12.or.us/secondary/howto/organize/images/timeline.gif
Here's one way:
# some fake data please, maestro
fakedata-sample(1:10,10)
# leave a bit of extra space beneath and to the left of the plot
par(mar=6,6,4,2)
# this function will probably end up in the plotrix package
time.line-function(x,y,at=NULL,labels=TRUE,tlticks=NULL,...) {
if(missing(x) missing (y))
stop(Usage: time.line(x,y,at=NULL,labels=TRUE))
plotrange-par(usr)
# get the graphical parameters
oldpar-par(no.readonly=TRUE)
# turn off clipping
par(xpd=TRUE)
if(missing(x)) {
# it's a horizontal line
segments(plotrange[1],y,plotrange[2],y,...)
ticklen-(plotrange[4]-plotrange[3])*0.02
if(!is.null(tlticks))
segments(tlticks,y+ticklen,tlticks,y-ticklen,...)
mwidth-strwidth(M)
# blank out the line where labels will appear
rect(at-mwidth,y-ticklen,at+mwidth,y+ticklen,col=white,border=FALSE)
# rotate the text
par(srt=90)
# draw the labels
text(at,y,labels,...)
}
if(missing(y)) {
# it's a vertical line
# draw the line
segments(x,plotrange[3],x,plotrange[4],...)
ticklen-(plotrange[2]-plotrange[1])*0.02
if(!is.null(tlticks))
segments(x+ticklen,tlticks,x-ticklen,tlticks,...)
mheight-strheight(M)
# blank out the line where labels will appear
rect(x-ticklen,at-mheight,x+ticklen,at+mheight,col=white,border=FALSE)
# draw the labels
text(x,at,labels,...)
}
# restore the parameters
par(oldpar)
}
# create a file with the positions and labels you want like this:
# 2.5,first
# 4,second
# 7,third
# 8.5,fourth
# call it labels.txt and read it in
tl.labels-read.table(labels.txt,sep=,)
plot(fakedata,xlab=)
# display a horizontal time line
time.line(x=-1,at=tl.labels[[1]],labels=tl.labels[[2]],col=black)
# now a vertical one
time.line(y=-1,at=tl.labels[[1]],labels=tl.labels[[2]],col=black)
__
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
__
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
Re: [R] embedding fonts in eps files
Rudi, If it turns out that fonts can not be embedded with R, then one option is to import/export the file through CorelDraw (or other vector drawing software equivalent). The 'export to eps' function in CorelDraw provides an option to embed all the fonts! It requires manual labour, but it will work. You will be able to embed any font available in Corel Draw: fonts galore!! Sander. Rudi Alberts wrote: Hi, I have to make eps files with fonts embedded. I use the following postscript command: postscript(fig3a.eps, width = 5.2756, height = 7.27, pointsize = 7,horizontal = FALSE, onefile = FALSE, paper = special,family = Times) plot(...) dev.off() Are fonts automatically embedded in this way? How can I see that? If not, how to do it? regards, Rudi. __ 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 __ 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
Re: [R] Time line plot in R? -- one more problem
Jim,
Inspired by the question about font embedding, I plotted the time line
script to a postscript file. To my disappointment, I can not make the
time line appear properly on the postscript graph. It seems that the
device does not know I have plotted something new below the original graph!?
Any suggesting how to resize the graph to plot the time line correctly
in postscript?
Thanks,
Sander.
Jim Lemon wrote:
Sander Oom wrote:
Dear R users,
In order to illustrate the possible effects of events on variables
plotted against time, I would like plot a time line of events along side
the plot of the variables.
The x-axis should be some time unit; the y-axis should be the variable
of interest; the time line should be plotted below the graph along the
same x-axis scale.
As I have many variables and different events lists, I would like to
write a script to read the events from a file and plot them together
with the other plot.
The time line should look something like this:
http://www.oslis.k12.or.us/secondary/howto/organize/images/timeline.gif
Here's one way:
# some fake data please, maestro
fakedata-sample(1:10,10)
# leave a bit of extra space beneath and to the left of the plot
par(mar=c(6,6,4,2))
# this function will probably end up in the plotrix package
time.line-function(x,y,at=NULL,labels=TRUE,tlticks=NULL,...) {
if(missing(x) missing (y))
stop(Usage: time.line(x,y,at=NULL,labels=TRUE))
plotrange-par(usr)
# get the graphical parameters
oldpar-par(no.readonly=TRUE)
# turn off clipping
par(xpd=TRUE)
if(missing(x)) {
# it's a horizontal line
segments(plotrange[1],y,plotrange[2],y,...)
ticklen-(plotrange[4]-plotrange[3])*0.02
if(!is.null(tlticks))
segments(tlticks,y+ticklen,tlticks,y-ticklen,...)
mwidth-strwidth(M)
# blank out the line where labels will appear
rect(at-mwidth,y-ticklen,at+mwidth,y+ticklen,col=white,border=FALSE)
# rotate the text
par(srt=90)
# draw the labels
text(at,y,labels,...)
}
if(missing(y)) {
# it's a vertical line
# draw the line
segments(x,plotrange[3],x,plotrange[4],...)
ticklen-(plotrange[2]-plotrange[1])*0.02
if(!is.null(tlticks))
segments(x+ticklen,tlticks,x-ticklen,tlticks,...)
mheight-strheight(M)
# blank out the line where labels will appear
rect(x-ticklen,at-mheight,x+ticklen,at+mheight,col=white,border=FALSE)
# draw the labels
text(x,at,labels,...)
}
# restore the parameters
par(oldpar)
}
# create a file with the positions and labels you want like this:
# 2.5,first
# 4,second
# 7,third
# 8.5,fourth
# call it labels.txt and read it in
tl.labels-read.table(labels.txt,sep=,)
plot(fakedata,xlab=)
# display a horizontal time line
time.line(x=-1,at=tl.labels[[1]],labels=tl.labels[[2]],col=black)
# now a vertical one
time.line(y=-1,at=tl.labels[[1]],labels=tl.labels[[2]],col=black)
__
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
__
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
RE: [R] IFELSE across large array?
Dear all,
Thanks to the help from the R mailing list we now have a script performing
the required tasks, i.e. applying a mask and filter to a 3D array. See
script below. Any further suggestions to simplify the code are very
welcome. We are currently debugging the larger script of which this section
is a small part.
The script could serve as a framework for many typical GIS neighborhood
analyses, in this case the 'block majority'. The 'var.range' variable is
derived from variography earlier in the main script.
Yours,
Sander and Alessandro
##INPUT SECTION1##
numRow-7
numCol-8
numReal-2
var.range- 2.1
cellsize- 0.25
##CALCULATE WINDOW SIZES BASED ON INPUT ABOVE#
#window size ...minimum is 1
numCells-round((var.range*cellsize)/2,0)
if(numCells 1) {numCells-1}
###
library(abind)
#Sim: 0=Absent; 1=Present; 10=NA
vectSim - c(1,1,0,0,1,1,10,1,0,10,1,1,1,0,1,1,1,1,1,0,0,1,0,0,10,0,0,
1,1,1,0,1,1,0,0,0,1,1,1,0,1,1,0,1,0,1,10,0,1,0,0,1,0,0,10,0,0,0,1,1,0,
0,1,0,1,0,0,1,0,0,1,1,1,1,1,10,1,10,1,0,1,0,1,0,0,1,1,0,10,10,1,1,0,1,
0,1,1,1,0,0,0,1,0,0,0,0,0,0,1,1,0,0)
#Mask: 10=NA; 1=DATA
vectMask - c(10,1,10,1,1,10,1,10,1,10,10,1,10,1,10,1,1,10,1,10,1,10,10,1,10,1,
10,1,1,10,1,10,1,10,10,1,10,1,10,1,1,10,1,10,1,10,10,1,10,1,10,1,1,10,1,10)
length(vectSim)
length(vectMask)
Sim - array(vectSim, c(numRow,numCol,numReal))
Sim[Sim==10]-NA
Sim
Mask - array(data=vectMask, c(numRow,numCol,numReal))
Mask
##expand array
junkcol-array(rep(NA,numReal*numRow),c(numRow,numCells,numReal))
#add columns borders
#cols
Sim-abind(junkcol,Sim,junkcol, along=2)
dim(Sim)[2]
junkrow-array(rep(NA,dim(Sim)[2]),c(numCells,dim(Sim)[2],numReal))
Sim-abind(junkrow,Sim,junkrow,along=1)
##DEFINE PORTION OF ARRAY ON WHICH MOVING WINDOW RUNS
##IT AVOIDS THE na BORDER
#maximum row and col indexes to consider are equal
# to num num rows and num cols in the original matrix
minr-1+numCells
maxr-dim(Sim)[1]-numCells
minc-1+numCells
maxc-dim(Sim)[2]-numCells
clean-function(a,nr,r,c) {
rmin-row(a)[r-numCells]
rmax-row(a)[r+numCells]
cmin-col(a)[nr*(c-numCells)]
cmax-col(a)[nr*(c+numCells)]
junk-a[rmin:rmax,cmin:cmax]
junk[numCells,numCells]-NA
junk-as.vector(junk)
sample(na.omit(junk),1)
}
for (n in 1:numReal) {
realiz-Sim[1:dim(Sim)[1],1:dim(Sim)[2],n]
rz-realiz
for (r in minr:maxr) {
for (c in minc:maxc) {
if(is.na(realiz[r,c]) ) {
realiz[r,c]-clean(a=realiz,nr=numRow,r=r,c=c)
Sim[1:dim(Sim)[1],1:dim(Sim)[2],n]-realiz
}#end if
}
}
}##end overall loop
# Return to original array dimensions, removing extra NA
minr - (numCells+1)
maxr - (dim(Sim)[1]-numCells)
minc - (numCells+1)
maxc - (dim(Sim)[2]-numCells)
expSim - Sim
expSim
Sim - expSim[minr:maxr, minc:maxc, 1:dim(Sim)[3]]
Sim
# Clip simulation results using the mask
Mask
Sim[Sim==10] - NA
Mask[Mask==10] - NA
A - Sim + Mask -1
A
#
At 01:55 2004/11/24, Ray Brownrigg wrote:
From: Liaw, Andy [EMAIL PROTECTED]
Date: Tue, 23 Nov 2004 12:28:48 -0500
I'll give it half a crack:
Steps a through c can be done via nested ifelse(), treating A and M as
vectors (as they really are). Step d is the hard one. I'd write a simple
Fortran code and use .Fortran() for that.
I don't see how any of the *apply() functions can help here, as your
operations are element-wise, not dimension-wise.
Andy
The original message mentions that the value 10 is actually NODATA,
and if one recodes 10 as NA, steps a) to c) become trivial, namely:
A[A == 10] - NA
M[M == 10] - NA
return(A + M - 1)
Then if step d) is performed first (i.e. appropriate values in A are
replaced by the 'most common neighbour' [perhaps using
round(mean(.., na.rm=T))] this still works, but would have to be repeated
for each replication (the third dimension).
Ray Brownrigg
From: Sander Oom
Dear all,
As our previous email did not get any response, we try again with a
reformulated question!
We are trying to do something which needs an efficient loop
over a huge
array, possibly functions such as apply and related (tapply,
lapply...?), but can't really understand syntax and examples in
practice...i.e. cant' make it work.
to be more specific:
we are trying to apply a mask to a 3D array.
By this I mean that when overlaying [i.e. comparing element
by element]
the mask on to the array the mask should change array
elements according to
the values of both array and mask elements
the mask has two values: 1 and 10.
the array elements have 3 values: 0, 1, or 10
sticking for the moment to the single 2d array case
for example:
[A= array] 100 10 1 10 0
1 10 1 0 0 10
[ M=mask] 1 10 10 1 1 1
101 1 1 10 10
I would like the array elements to:
a) IF A(ij) !=10 and Mij = 1
leave A(ij
[R] IFELSE across large array?
Dear all, As our previous email did not get any response, we try again with a reformulated question! We are trying to do something which needs an efficient loop over a huge array, possibly functions such as apply and related (tapply, lapply...?), but can't really understand syntax and examples in practice...i.e. cant' make it work. to be more specific: we are trying to apply a mask to a 3D array. By this I mean that when overlaying [i.e. comparing element by element] the mask on to the array the mask should change array elements according to the values of both array and mask elements the mask has two values: 1 and 10. the array elements have 3 values: 0, 1, or 10 sticking for the moment to the single 2d array case for example: [A= array] 100 10 1 10 0 1 10 1 0 0 10 [ M=mask] 1 10 10 1 1 1 101 1 1 10 10 I would like the array elements to: a) IF A(ij) !=10 and Mij = 1 leave A(ij) unchanged b) IF A(ij) != 10 and M(ij) =10 change A(ij) to M(ij) i.e mask value (10) c)IF A(ij) = 10 and M(ij) = 10 leave (Aij) unchanged d) IF A(ij) = 10 and M(ij) !=10 replace A(ij) with the majority value in the 8-neighborhood (or whatever if it is an edge element) BUT ignoring 10s in this neighborhood (i.e. with either 1 or 0, whichever is in majority) because the array is 3d I would like to repeat the thing with all the k elements (2d sub-arrays) of the array in question, using the same mask for al k elements Would you be able to suggest a strategy to do this? thanks very much Alessandro and Sander. __ [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
[R] IFELSE across a 3D array?
Dear all, We are trying to clean multiple realizations of a pattern. Erroneous NODATA and spurious DATA occur in the realizations. As we have to do a 1000 realizations for many patterns, efficiency of the code is important. We need to correct the realizations with a 'mask' pattern of DATA/NODATA. We think an ifelse should do the job. Spurious DATA will be simply removed using the mask. Erroneous NODATA should be filled in by averaging across the cell's nearest neighbors. The ifelse table is as follows: IFELSE criteria to determine TARGET from realization and mask Real Mask Target 10 10 10 10 1 0/1 determined through post processing function (average across neighbors) 0/110 10 0/1 1 0/1 We think that an APPLY on a multi dimensional array is the way to go, with each realization (2 dimensions) being a dimension of the array (like a stack of maps) This is where we have got so far: * library(abind) #Sim: 0=Absent; 1=Present; 10=NODATA vectSim - c(0,1,0,1,10,0,1,0,1,1,10,0,1,1,0,1,0,10,10,1,0,1,0,1,0,1,1,10,0,10,10,1,0,1,0,10) #Mask: 1=DATA; 10=NODATA vectMask - c(10,1,10,1,1,10,1,10,1,10,10,1) length(vectSim) length(vectMask) numRow-3 numCol-4 numReal-3 Sim - array(vectSim, c(numRow,numCol,numReal)) Sim Mask - array(vectMask, c(numRow,numCol)) Mask SmoothSim - apply(Sim, c(1,2), ifelse(MaskSim==10,33,99)) SmoothSim Any help is much appreciated! Thanks, Sander and Alessandro. __ [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
[R] Calculating distances between points in a data frame?
Dear list,
I would like to calculate the distance between consecutive points in a data
frame. Of course the first point in the data frame does not have a point of
origin, and should get a value NA. I have tried two different loops, which
both result in error:
num - seq(0,10,1)
X - seq(0,30,3)
Y - seq(0,40,4)
XY - data.frame(num, X, Y)
attach(XY)
summary(XY)
num X Y
Min. : 0.0 Min. : 0.0 Min. : 0
1st Qu.: 2.5 1st Qu.: 7.5 1st Qu.:10
Median : 5.0 Median :15.0 Median :20
Mean : 5.0 Mean :15.0 Mean :20
3rd Qu.: 7.5 3rd Qu.:22.5 3rd Qu.:30
Max. :10.0 Max. :30.0 Max. :40
plot(X,Y)
rngNum - range(num)
for (i in rngNum){
+ XY$DistXY[i] - sqrt( ((X[i]-X[i-1])^2) + ((Y[i]-Y[i-1])^2) )
+ }
Error in $-.data.frame(`*tmp*`, DistXY, value = sqrt(((X[i] - X[i - :
replacement has 10 rows, data has 11
for (i in rngNum){
+ XY$DistXY2[i] - ifelse(i=min(rngNum), NA, sqrt(((X[i]-X[i-1])^2) +
((Y[i]-Y[i-1])^2)) )
+ }
Error in ifelse(i = min(rngNum), NA, sqrt(((X[i] - X[i - 1])^2) + ((Y[i] - :
unused argument(s) (i ...)
detach(XY)
Any suggestions much appreciated,
Sander Oom.
--
Dr. Sander P. Oom
Animal, Plant and Environmental Sciences
University of the Witwatersrand
Private Bag 3
Wits 2050
South Africa
Tel (work) +27 (0)11 717 64 04
Tel (home) +27 (0)18 297 44 51
Fax +27 (0)18 299 24 64
Email [EMAIL PROTECTED]
Web www.oomvanlieshout.net/sander
__
[EMAIL PROTECTED] mailing list
https://www.stat.math.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] Calculating distances between points in a data frame?
Hi Gabor,
Thanks for your suggestion. However when installing the package gregmisc, I
get the following error:
local({a - CRAN.packages()
+ install.packages(select.list(a[,1],,TRUE), .libPaths()[1], available=a)})
trying URL `http://cran.r-project.org/bin/windows/contrib/1.9/PACKAGES'
Content type `text/plain; charset=iso-8859-1' length 17940 bytes
opened URL
downloaded 17Kb
trying URL
`http://cran.r-project.org/bin/windows/contrib/1.9/gregmisc_1.11.0.zip'
Error in download.file(url, destfile, method, mode = wb) :
cannot open URL
`http://cran.r-project.org/bin/windows/contrib/1.9/gregmisc_1.11.0.zip'
In addition: Warning message:
cannot open: HTTP status was `404 Not Found'
I was quite surprised as well! I tried different mirrors, but to no avail.
Maybe tomorrow,
Sander.
At 18:44 2004/05/31, you wrote:
Try using running from the gregmisc package with pad = TRUE:
require(gregmisc)
XY - data.frame(num = seq(0,10), X = seq(0,30,3), Y = seq(0, 40, 4) )
DistXY - function(idx) {
i - idx[2]
with(XY, sqrt( (X[i]-X[i-1])^2 + (Y[i]-Y[i-1])^2 ) )
}
XY$Dist - running( 1:nrow(XY), width=2, fun = DistXY, pad = TRUE )
Sander Oom slist at oomvanlieshout.net writes:
:
: Dear list,
:
: I would like to calculate the distance between consecutive points in a data
: frame. Of course the first point in the data frame does not have a point of
: origin, and should get a value NA. I have tried two different loops, which
: both result in error:
:
: num - seq(0,10,1)
: X - seq(0,30,3)
: Y - seq(0,40,4)
: XY - data.frame(num, X, Y)
: attach(XY)
: summary(XY)
:num X Y
: Min. : 0.0 Min. : 0.0 Min. : 0
: 1st Qu.: 2.5 1st Qu.: 7.5 1st Qu.:10
: Median : 5.0 Median :15.0 Median :20
: Mean : 5.0 Mean :15.0 Mean :20
: 3rd Qu.: 7.5 3rd Qu.:22.5 3rd Qu.:30
: Max. :10.0 Max. :30.0 Max. :40
: plot(X,Y)
: rngNum - range(num)
: for (i in rngNum){
: + XY$DistXY[i] - sqrt( ((X[i]-X[i-1])^2) + ((Y[i]-Y[i-1])^2) )
: + }
: Error in $-.data.frame(`*tmp*`, DistXY, value = sqrt(((X[i] - X[i - :
: replacement has 10 rows, data has 11
: for (i in rngNum){
: + XY$DistXY2[i] - ifelse(i=min(rngNum), NA, sqrt(((X[i]-X[i-1])^2) +
: ((Y[i]-Y[i-1])^2)) )
: + }
: Error in ifelse(i = min(rngNum), NA, sqrt(((X[i] - X[i - 1])^2) + ((Y[i]
- :
: unused argument(s) (i ...)
: detach(XY)
:
:
: Any suggestions much appreciated,
:
: Sander Oom.
:
: --
: Dr. Sander P. Oom
: Animal, Plant and Environmental Sciences
: University of the Witwatersrand
: Private Bag 3
: Wits 2050
: South Africa
:
: Tel (work) +27 (0)11 717 64 04
: Tel (home) +27 (0)18 297 44 51
: Fax +27 (0)18 299 24 64
:
: Email sander at oomvanlieshout.net
: Web www.oomvanlieshout.net/sander
:
: __
: R-help at stat.math.ethz.ch mailing list
: https://www.stat.math.ethz.ch/mailman/listinfo/r-help
: PLEASE do read the posting guide!
http://www.R-project.org/posting-guide.html
:
:
__
[EMAIL PROTECTED] mailing list
https://www.stat.math.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
__
[EMAIL PROTECTED] mailing list
https://www.stat.math.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
RE: Re: [R] Plotting Time against Date for time series data?
The plot is nearly there! Using the axis.POSIXct command I have got the
x-axis under control. However, the units for the y-axis (Time) are in
seconds by default (i.e. range is from 0 to 1440). I'm trying to plot hours
along the y-axis, without changing the units for the plot itself, but
without any luck.
The #'s were a quick solution to stop WinEdt reformatting my code
automatically, now I discovered I can switch to a different wrap mode.
This is the code just now:
psDateTime - function(DateTime, FileName){
+ strFileName - paste(Analysis\\Graphics\\time_date_, FileName, sep
= )
+ startMonth - as.POSIXct(as.Date(01/09/2003, format=%d/%m/%Y))
+ endMonth - as.POSIXct(as.Date(01/06/2004, format=%d/%m/%Y))
+ startTime - as.POSIXct(as.Date(00:00:00, format=%H:%m:%s))
+ endTime - as.POSIXct(as.Date(23:59:59, format=%H:%m:%s))
+ Date - trunc(DateTime, day)
+ Time - DateTime - Date
+ postscript(strFileName)
+ plot(x=Date, y=Time,
+ xlab= Date (month/year), ylab= Time (hours), xaxt=n, #yaxt=n,
+ xlim=c(as.POSIXct(startMonth), as.POSIXct(endMonth)),
+ ylim=c(as.POSIXct(startTime), as.POSIXct(endTime))
+ )
+ axis.POSIXct(1, at=seq(startMonth, endMonth, by=month), format=%m/%y)
+ axis.POSIXct(2, at=seq(startTime, endTime, by=hours), format=%H)
+ pstamp(pwd=FALSE, time=TRUE)
+ dev.off()
+ }
psDateTime(datKruger1$DaTim, datKruger1.eps)
Error in plot.window(xlim, ylim, log, asp, ...) :
need finite ylim values
and when removing the ylim line from the above code, an error associated
with the 'axis.POSICct(2)' line:
psDateTime(datKruger1$DaTim, datKruger1.eps)
Error in if (to = from) stop(`to' must be later than `from') :
missing value where TRUE/FALSE needed
Any help much appreciated,
Sander.
At 15:58 2004/05/17, you wrote:
On Mon, 17 May 2004, Slist wrote:
I have a data set containing GPS fixes of animal locations. To check that
the GPS's are working properly, I would like to plot the time of the fixes
(y-axis) against the date of the fixes (x-axis). If all works well, the
plot should show four regular fixes per day. The x-axis should be labelled
with month/year (i.e. 11/04) and the y-axis by hour from 00 to 24. I would
like to control the x-axis limits.
I have looked at several date and time related help pages, but get
horribly
lost in all the terminology. The main challenge is to isolate date and
time
from the date/time object for plotting (marked ???). Therefore, I would
like the following example code to work:
dates - c(02/27/92, 02/27/92, 01/14/92, 01/14/92, 03/28/92,
03/28/92) #
times - c(23:03:20, 10:29:56, 01:03:30, 13:03:30, 18:21:03,
06:56:26) #
x - paste(dates, times) #
DateTime - strptime(x, %m/%d/%y %H:%M:%S) #
Why do you end lines with #? It is rather confusing.
Date - trunc(DateTime, day)
Time - DateTime - Date
plot(Date, Time)
appears to do what you want (except the x axis labelling, which you can
alter by a call to axis.POSIXct and x-axis limits, which need to be set
via xlim as a POSIXct object.).
--
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
__
[EMAIL PROTECTED] mailing list
https://www.stat.math.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
[R] Plotting Time against Date for time series data? (2)
Almost but not quite:
Someone told me R would produce such high quality graphics, I would never
need a separate graphics package again. This does however mean that I need
to be able to draw the graph exactly the way I want it to look!
I have searched the web and help files extensively, but there are few
references on plotting time on an axis. The function presented here
(http://tolstoy.newcastle.edu.au/R/help/00a/1031.html) could be useful, but
I'm not sure how to implement it in my example. Anyway there should be a
simpler way.
The following code uses the example again. It does have hours on the
x-axis, but the tic marks are printed as: 0, 5, 10, 15, 20. Not a very
intuitive frequency of tic marks. How do I get R to plot tic marks at: 0,
4, 8, 12, 16, 20, 24 or 0, 2, 4..., 24? I want the tic marks to be fixed,
i.e. they should still go from 0-24 even when the range of the data is smaller.
dates - c(02/27/92, 02/27/92, 01/14/92, 01/14/92, 03/28/92,
03/28/92) #
times - c(23:03:20, 10:29:56, 01:03:30, 13:03:30, 18:21:03,
06:56:26) #
x - paste(dates, times) #
DateTime - strptime(x, %m/%d/%y %H:%M:%S) #
Date - trunc(DateTime, day)
Time - difftime(DateTime, Date, units=hours)
startMonth - as.POSIXct(as.Date(01/01/1992, format=%d/%m/%Y))
endMonth - as.POSIXct(as.Date(01/05/1992, format=%d/%m/%Y))
plot(x=Date, y=Time,
xlab= Date (month/year), ylab= Time (hours), xaxt=n, yaxt=n,
xlim=c(startMonth, endMonth),
ylim=c(0, 24)
) #
axis.POSIXct(1, at=seq(startMonth, endMonth, by=month), format=%m/%y)
#r - as.POSIXct(range(Time), hours) #
clock24 -
strptime(c(0,2,4,6,8,10,12,14,16,18,20,22,23), %H)
r - as.POSIXct(round(range(clock24), hours))
axis.POSIXct(2, at=seq(r[1], r[2], by=hour), format=%H) #
More suggestions welcome,
Sander
At 18:39 2004/05/17, you wrote:
On Mon, 17 May 2004, Sander Oom wrote:
The plot is nearly there! Using the axis.POSIXct command I have got the
x-axis under control. However, the units for the y-axis (Time) are in
seconds by default (i.e. range is from 0 to 1440). I'm trying to plot
hours
along the y-axis, without changing the units for the plot itself, but
without any luck.
I got hours for your example, I believe, but it does depend on the data.
Try
Time - difftime(Date, DateTime, units=hours)
and treat them as numbers, not as POSIXct (they are not dates).
The #'s were a quick solution to stop WinEdt reformatting my code
automatically, now I discovered I can switch to a different wrap mode.
This is the code just now:
psDateTime - function(DateTime, FileName){
+ strFileName - paste(Analysis\\Graphics\\time_date_, FileName, sep
= )
+ startMonth - as.POSIXct(as.Date(01/09/2003, format=%d/%m/%Y))
+ endMonth - as.POSIXct(as.Date(01/06/2004, format=%d/%m/%Y))
+ startTime - as.POSIXct(as.Date(00:00:00, format=%H:%m:%s))
+ endTime - as.POSIXct(as.Date(23:59:59, format=%H:%m:%s))
Just startTime - 0, endTime - 24 should do.
+ Date - trunc(DateTime, day)
+ Time - DateTime - Date
+ postscript(strFileName)
+ plot(x=Date, y=Time,
+ xlab= Date (month/year), ylab= Time (hours), xaxt=n,
#yaxt=n,
+ xlim=c(as.POSIXct(startMonth), as.POSIXct(endMonth)),
+ ylim=c(as.POSIXct(startTime), as.POSIXct(endTime))
+ )
+ axis.POSIXct(1, at=seq(startMonth, endMonth, by=month),
format=%m/%y)
+ axis.POSIXct(2, at=seq(startTime, endTime, by=hours), format=%H)
+ pstamp(pwd=FALSE, time=TRUE)
+ dev.off()
+ }
psDateTime(datKruger1$DaTim, datKruger1.eps)
Error in plot.window(xlim, ylim, log, asp, ...) :
need finite ylim values
and when removing the ylim line from the above code, an error associated
with the 'axis.POSICct(2)' line:
psDateTime(datKruger1$DaTim, datKruger1.eps)
Error in if (to = from) stop(`to' must be later than `from') :
missing value where TRUE/FALSE needed
Any help much appreciated,
Sander.
At 15:58 2004/05/17, you wrote:
On Mon, 17 May 2004, Slist wrote:
I have a data set containing GPS fixes of animal locations. To
check that
the GPS's are working properly, I would like to plot the time of
the fixes
(y-axis) against the date of the fixes (x-axis). If all works well, the
plot should show four regular fixes per day. The x-axis should be
labelled
with month/year (i.e. 11/04) and the y-axis by hour from 00 to 24.
I would
like to control the x-axis limits.
I have looked at several date and time related help pages, but get
horribly
lost in all the terminology. The main challenge is to isolate date and
time
from the date/time object for plotting (marked ???). Therefore, I would
like the following example code to work:
dates - c(02/27/92, 02/27/92, 01/14/92, 01/14/92, 03/28/92,
03/28/92) #
times - c(23:03:20, 10:29:56, 01:03:30, 13:03:30, 18:21:03,
06:56:26) #
x - paste(dates, times) #
DateTime - strptime(x, %m/%d/%y %H:%M:%S
