Re: [R] Aggregate rainfall data

2016-07-19 Thread roslinazairimah zakaria
Hi David,

Thank you so much for your help and others.  Here is the code.

balok <- read.csv("G:/A_backup 11 mei 2015/DATA (D)/1 Universiti Malaysia
Pahang/ISM-3 2016 UM/Data/Hourly Rainfall/balok2.csv",header=TRUE)
head(balok, 10); tail(balok, 10)
str(balok)

## Introduce NAs for
balok$Rain.mm2 <- as.numeric( as.character(balok$Rain.mm) )
head(balok$Rain.mm2); tail(balok$Rain.mm2)
head(balok, 10); tail(balok, 10)

## Change date format from DD/MM/ to Day, Month, Year separately
realdate <- as.Date(balok$Date,format="%d/%m/%Y")
dfdate <- data.frame(date=realdate)
year=as.numeric (format(realdate,"%Y"))
month=as.numeric (format(realdate,"%m"))
day=as.numeric (format(realdate,"%d"))

balok2 <-cbind(dfdate,day,month,year,balok)
colnames(balok2)
head(balok2)

## New data format
balok2_new <- balok2[,c(-1,-5,-7)]
colnames(balok2_new)
head(balok2_new); tail(balok2_new)

## Aggregate data
## Sum rainfall amount from HOURLY to DAILY
dt <- balok2_new
str(dt)
aggbalok_day <- aggregate(dt[,5], by=dt[,c(1,2,3)],FUN=sum, na.rm=TRUE)
head(aggbalok_day)

## Sum rainfall amount from HOURLY to MONTHLY
dt <- balok2_new
str(dt)
aggbalok_mth <- aggregate(dt[,5], by=dt[,c(2,3)],FUN=sum, na.rm=TRUE)
head(aggbalok_mth)

Now I would like to find the basic statistics summary for the data
according to monthly.


Best regards




On Wed, Jul 13, 2016 at 10:37 PM, David Winsemius 
wrote:

>
> > On Jul 13, 2016, at 3:21 AM, roslinazairimah zakaria <
> roslina...@gmail.com> wrote:
> >
> > Dear David,
> >
> > I got your point.  How do I remove the data that contain "0.0?".
> >
> > I tried : balok <- cbind(balok3[,-5],
> balok3$Rain.mm[balok3$Rain.mm==0.0?] <- NA)
>
> If you had done as I suggested, the items with factor levels of "0.0?"
> would have automatically become NA and you would have gotten a warning
> message:
>
> > testfac <- factor( c(rep("0.0",7), "0.07", "0.0?", '0.01', '0.17'))
> > testfac
>  [1] 0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.07 0.0? 0.01 0.17
> Levels: 0.0 0.0? 0.01 0.07 0.17
> > as.numeric(as.character( testfac))
>  [1] 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.07   NA 0.01 0.17
> Warning message:
> NAs introduced by coercion
>
>
>
> >
> > However all the Rain.mm column all become NA.
> >
> >day month year Time balok3$Rain.mm[balok3$Rain.mm == "0.0?"] <- NA
> > 1   30 7 2008  9:00:00 NA
> > 2   30 7 2008 10:00:00 NA
> > 3   30 7 2008 11:00:00 NA
> > 4   30 7 2008 12:00:00 NA
> > 5   30 7 2008 13:00:00 NA
> > 6   30 7 2008 14:00:00 NA
> > 7   30 7 2008 15:00:00 NA
> > 8   30 7 2008 16:00:00 NA
> > 9   30 7 2008 17:00:00 NA
> > 10  30 7 2008 18:00:00 NA
> >
> > Thank you so much.
> >
> >
> > On Wed, Jul 13, 2016 at 9:42 AM, David Winsemius 
> wrote:
> >
> > > On Jul 12, 2016, at 3:45 PM, roslinazairimah zakaria <
> roslina...@gmail.com> wrote:
> > >
> > > Dear R-users,
> > >
> > > I have these data:
> > >
> > > head(balok, 10); tail(balok, 10)
> > >Date Time Rain.mm
> > > 1  30/7/2008  9:00:00   0
> > > 2  30/7/2008 10:00:00   0
> > > 3  30/7/2008 11:00:00   0
> > > 4  30/7/2008 12:00:00   0
> > > 5  30/7/2008 13:00:00   0
> > > 6  30/7/2008 14:00:00   0
> > > 7  30/7/2008 15:00:00   0
> > > 8  30/7/2008 16:00:00   0
> > > 9  30/7/2008 17:00:00   0
> > > 10 30/7/2008 18:00:00   0
> > >   Date Time Rain.mm
> > > 63667 4/11/2015  3:00:00   0
> > > 63668 4/11/2015  4:00:00   0
> > > 63669 4/11/2015  5:00:00   0
> > > 63670 4/11/2015  6:00:00   0
> > > 63671 4/11/2015  7:00:00   0
> > > 63672 4/11/2015  8:00:00   0
> > > 63673 4/11/2015  9:00:00 0.1
> > > 63674 4/11/2015 10:00:00 0.1
> > > 63675 4/11/2015 11:00:00 0.1
> > > 63676 4/11/2015 12:00:000.1?
> > >
> > >> str(balok)
> > > 'data.frame':   63676 obs. of  3 variables:
> > > $ Date   : Factor w/ 2654 levels "1/1/2009","1/1/2010",..: 2056 2056
> 2056
> > > 2056 2056 2056 2056 2056 2056 2056 ...
> > > $ Time   : Factor w/ 24 levels "1:00:00","10:00:00",..: 24 2 3 4 5 6 7
> 8 9
> > > 10 ...
> > > $ Rain.mm: Factor w/ 352 levels "0","0.0?","0.1",..: 1 1 1 1 1 1 1 1 1
> 1
> >
> > Thar's your problem:
> >
> >   Rain.mm: Factor w/ 352 levels "0","0.0?","0.1"
> >
> > Need to use the standard fix for the screwed-up-factor-on-input-problem
> >
> >   balok$Rain.mm2 <- as.numeric( as.character(balok$Rain.mm) )
> >
> > Cannot just do as.numeric because factors are actually already numeric.
> >
> > --
> > David.
> >
> >
> > > ...
> > >
> > > and I 

Re: [R] Build command in library(devtools)

2016-07-19 Thread Steven Yen
Thanks. I found the reason was Rtools does not run under the new version 
of R. I had to go back to as early as R 3.0.2 (September 2013) to make 
Rtools work.
Any idea for a go-around? Thanks.

On 7/19/2016 4:38 PM, John McKown wrote:
> On Tue, Jul 19, 2016 at 3:15 PM, Steven Yen  >wrote:
>
> I recently updated my R and RStudio to the latest version and now the
> binary option in the "build" command in devtools stops working.
>
> I went around and used the binary=F option which worked by I get the
> .tar.gz file instead of the .zip file which I prefer.
>
> Does anyone understand the following error message:
>
> status 127
> running 'zip' failed
>
>
> ​I'm not totally sure, but I think that means that R cannot find the 
> "zip" program in order to run it. ​
>
>
> ===
> setwd("A:/R/yenlib/"); library(devtools)
> #build("yenlib",binary=T) # Thisfailed with an error message
> build("yenlib",binary=F) # This works
>
>  > build("yenlib",binary=T)
> "C:/PROGRA~1/R/R-33~1.1/bin/x64/R" --no-site-file  \
>--no-environ --no-save --no-restore --quiet CMD INSTALL \
>"A:\R\yenlib\yenlib" --build
>
> * installing to library
> 'C:/Users/syen01/AppData/Local/Temp/Rtmp8A7KEw/temp_libpath4074149a528e'
> * installing *source* package 'yenlib' ...
> ** R
> ** data
> ** preparing package for lazy loading
> ** help
> *** installing help indices
> ** building package indices
> ** testing if installed package can be loaded
> *** arch - i386
> *** arch - x64
> * MD5 sums
> Warning: running command '"zip" -r9Xq "A:/R/yenlib/yenlib_16.3.zip"
> yenlib' had status 127
> running 'zip' failed
> * DONE (yenlib)
> [1] "A:/R/yenlib/yenlib_16.3.zip"
>  >
>
>
> [[alternative HTML version deleted]]
>
> __
> R-help@r-project.org  mailing list --
> To UNSUBSCRIBE and more, see
> 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.
>
>
>
>
> -- 
> "Worry was nothing more than paying interest on a loan that a man may 
> never borrow"
>
> From: "Quest for the White Wind" by Alan Black
>
> Maranatha! <><
> John McKown


[[alternative HTML version deleted]]

__
R-help@r-project.org mailing list -- To UNSUBSCRIBE and more, see
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] pairs: adjusting margins and labeling axes

2016-07-19 Thread Duncan Mackay
Hi
Will doing in lattice suite

>From  https://stat.ethz.ch/pipermail/r-help/2007-October/142116.html and
https://stat.ethz.ch/pipermail/r-help/2007-October/142124.html]
This is a direct cut and paste from the last url to give you an idea of
Deepayan Sarkar's script

library(lattice)

panel.corval2 <- function(x, y, digits=2, prefix="", cex.cor, ...)
{
require(grid)
r <- abs(cor(x, y, use = "complete.obs"))
txt <- format(c(r, 0.123456789), digits=digits)[1]
txt <- paste(prefix, txt, sep="")
if (missing(cex.cor)) cex.cor <- 10 / nchar(txt)
grid.text(txt, 0.5, 0.5, gp = gpar(cex = cex.cor))
}

splom(iris[1:4], groups = iris$Species, pch = 16,
  lower.panel = function(...) {
  panel.xyplot(...)
  panel.loess(..., col = 1, lwd = 3)
  },
  upper.panel = panel.corval2)


Regards

Duncan

Duncan Mackay
Department of Agronomy and Soil Science
University of New England
Armidale NSW 2351
Email: home: mac...@northnet.com.au

-Original Message-
From: R-help [mailto:r-help-boun...@r-project.org] On Behalf Of michael
young
Sent: Wednesday, 20 July 2016 03:30
To: r-help@r-project.org
Subject: [R] pairs: adjusting margins and labeling axes

The default shape for this correlation scatterplot is rectangle.  I changed
it to square, but then the x-axis spacing between squares are off.  Is
there an easy way to change x-axis spacing between squares to that of the
y-axis spacing size?

I decided to hide the name values of the diagonal squares.  I want them
along the x and y axis instead, outside of the fixed number scale I have.
I haven't seen any online example of 'pairs' with this and all my searches
have yielded nothing.  Any ideas?  Thanks

par(pty="s")
panel.cor <- function(x, y, digits = 2, prefix="", cex.cor, ...)
{
usr <- par("usr"); on.exit(par(usr))
par(usr = c(0, 1, 0, 1),xlog=FALSE,ylog=FALSE)
# correlation coefficient
r <- cor(x, y)
txt <- format(c(r, 0.123456789), digits = digits)[1]
txt <- paste("r= ", txt, sep = "")
if(missing(cex.cor)) cex.cor <- 0.8/strwidth(txt)
text(0.5, 0.6, txt, cex=cex.cor * r)

# p-value calculation
p <- cor.test(x, y)$p.value
txt2 <- format(c(p, 0.123456789), digits = digits)[1]
txt2 <- paste("p= ", txt2, sep = "")
if(p<0.01) txt2 <- paste("p= ", "<0.01", sep = "")
text(0.5, 0.4, txt2)
}

pairs(iris, upper.panel = panel.cor,xlim=c(0.1,10),
ylim=c(0.1,10),log="xy",text.panel = NULL,pch=".")

[[alternative HTML version deleted]]

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

__
R-help@r-project.org mailing list -- To UNSUBSCRIBE and more, see
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] Sampe numbers

2016-07-19 Thread Daniel Nordlund

Frederic,

you need to be more clear about what you want to do.  If you are 
sampling 100 numbers from 1:100 with replacement, the sum will always be 
greater than 50 (as has already been pointed out).


In addition, you mention trying to eliminate zeros.  If you are sampling 
from 1:100 there will be no zeros.  So what are you really trying to do?



Dan

Dan Nordlund
Port Townsend, WA  USA


On 7/19/2016 11:54 AM, Frederic Ntirenganya wrote:

The problem is that there
are some zeros  while the numbers should range from 1 to  100

As the replace= TRUE if the smallest candidate is 1 may be replaced
because if I can have 1, fifty times this is equal to 50.
Remember that it is a sampling with a replacement

I am still wondering how I can remove the zero value.


Frederic Ntirenganya
Maseno University,
African Maths Initiative,
Kenya.
Mobile:(+254)718492836
Email: fr...@aims.ac.za
https://sites.google.com/a/aims.ac.za/fredo/

On Tue, Jul 19, 2016 at 9:11 PM,  wrote:


N <- 100
C <- 50
x <- numeric(N)
for (i in 1:N){
  x[i] <- sample(C-sum(x),1)
}
x
sum(x)




*Frederic Ntirenganya >*
Sent by: "R-help" 

07/19/2016 01:41 PM
To
"r-help@r-project.org" ,

cc
Subject
[R] Sampe numbers




Hi Guys,
I am trying to sample 100 numbers from 1:100
i did it  like this:

sample(1:100,100, replace= TRUE)
but i want again include a constraint that their sum must be equal to 50

How could I do it?

Cheers

Frederic Ntirenganya
Maseno University,
African Maths Initiative,
Kenya.
Mobile:(+254)718492836
Email: fr...@aims.ac.za
https://sites.google.com/a/aims.ac.za/fredo/


[[alternative HTML version deleted]]

__
R-help@r-project.org mailing list -- To UNSUBSCRIBE and more, see
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide
http://www.R-project.org/posting-guide.html

and provide commented, minimal, self-contained, reproducible code.




[[alternative HTML version deleted]]

__
R-help@r-project.org mailing list -- To UNSUBSCRIBE and more, see
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.




--
Daniel Noredlund
Bothell, WA USA

__
R-help@r-project.org mailing list -- To UNSUBSCRIBE and more, see
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] Build command in library(devtools)

2016-07-19 Thread John McKown
On Tue, Jul 19, 2016 at 3:15 PM, Steven Yen  wrote:

> I recently updated my R and RStudio to the latest version and now the
> binary option in the "build" command in devtools stops working.
>
> I went around and used the binary=F option which worked by I get the
> .tar.gz file instead of the .zip file which I prefer.
>
> Does anyone understand the following error message:
>
> status 127
> running 'zip' failed
>

​I'm not totally sure, but I think that means that R cannot find the "zip"
program in order to run it. ​



>
> ===
> setwd("A:/R/yenlib/"); library(devtools)
> #build("yenlib",binary=T) # Thisfailed with an error message
> build("yenlib",binary=F) # This works
>
>  > build("yenlib",binary=T)
> "C:/PROGRA~1/R/R-33~1.1/bin/x64/R" --no-site-file  \
>--no-environ --no-save --no-restore --quiet CMD INSTALL \
>"A:\R\yenlib\yenlib" --build
>
> * installing to library
> 'C:/Users/syen01/AppData/Local/Temp/Rtmp8A7KEw/temp_libpath4074149a528e'
> * installing *source* package 'yenlib' ...
> ** R
> ** data
> ** preparing package for lazy loading
> ** help
> *** installing help indices
> ** building package indices
> ** testing if installed package can be loaded
> *** arch - i386
> *** arch - x64
> * MD5 sums
> Warning: running command '"zip" -r9Xq "A:/R/yenlib/yenlib_16.3.zip"
> yenlib' had status 127
> running 'zip' failed
> * DONE (yenlib)
> [1] "A:/R/yenlib/yenlib_16.3.zip"
>  >
>
>
> [[alternative HTML version deleted]]
>
> __
> R-help@r-project.org mailing list -- To UNSUBSCRIBE and more, see
> 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.
>



-- 
"Worry was nothing more than paying interest on a loan that a man may never
borrow"

From: "Quest for the White Wind" by Alan Black

Maranatha! <><
John McKown

[[alternative HTML version deleted]]

__
R-help@r-project.org mailing list -- To UNSUBSCRIBE and more, see
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] Label multiple-line plot

2016-07-19 Thread John Wasige
Thanks so much
Best Rgds
John

On Tue, Jul 19, 2016 at 9:56 PM, Greg Snow <538...@gmail.com> wrote:

> with 72 lines to label it will be crowded whatever you do.
>
> Here is one option (though I am showing with fewer lines):
>
> x <- sapply(1:15, function(i) cumsum(rnorm(100)))
> par(mar=c(5,4,4,3)+0.1)
> matplot(1:100, x, type='l', lty=1)
> mtext(LETTERS[1:15], side=4, at=x[100,], las=1, line=1)
>
> One way to spread out the labels is:
>
> library(TeachingDemos)
> new.y <- spread.labs(x[100,], mindiff = 1.2*strheight('A'),
>  min=min(x), max=max(x))
> mtext(LETTERS[1:15], side=4, at=new.y, las=1, line=2)
>
>
> There is also another label spreading function in the plotrix package
> and another option for finding space in a plot for labeling lines in
> either the rms or Hmisc package.  There are probably other packages
> with tools available now as well, but those are the ones that I am
> familiar with.
>
>
>
> On Tue, Jul 19, 2016 at 1:36 PM, John Wasige  wrote:
> > Thanks! It worked  for me.
> > matplot(for_jhon$ID, for_jhon[,2:73], type='l')
> >
> > Any I dea on how I can label multiple-line plot based on column names?
> >
> > Thanks for your help
> >
> > John
> >
> > On Tue, Jul 19, 2016 at 8:46 PM, Greg Snow <538...@gmail.com> wrote:
> >>
> >> Most attachments get stripped off, so your data did not make it through.
> >>
> >> But try:
> >>
> >> matplot(for_jhon$ID, for_jhon[,2:73], type='l')
> >>
> >>
> >> On Tue, Jul 19, 2016 at 12:24 PM, John Wasige 
> >> wrote:
> >> > Dear all,
> >> >
> >> > This is to kindly request for your help. I would like to plot my data.
> >> >
> >> > The R script below gives some plot that is not clear. How can I get a
> >> > clear
> >> > multiple-line plot. The data is attached herewith.
> >> >
> >> > ##R Script
> >> > for_jhon = read.csv("C:/LVM_share/for_ jhon.csv", header=TRUE,
> sep=";")
> >> > matplot(for_jhon$ID, cbind(for_jhon[,2:73]))
> >> >
> >> > Thanks for your help
> >> >
> >> > John
> >> > __
> >> > R-help@r-project.org mailing list -- To UNSUBSCRIBE and more, see
> >> > 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.
> >>
> >>
> >>
> >> --
> >> Gregory (Greg) L. Snow Ph.D.
> >> 538...@gmail.com
> >
> >
> >
> >
> > --
> > John Wasige
> > "There are no REGRATES in LIFE, just lessons (Jennifer Aniston)”
>
>
>
> --
> Gregory (Greg) L. Snow Ph.D.
> 538...@gmail.com
>



-- 
John Wasige
"There are no REGRATES in LIFE, just lessons (Jennifer Aniston)”

[[alternative HTML version deleted]]

__
R-help@r-project.org mailing list -- To UNSUBSCRIBE and more, see
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] Build command in library(devtools)

2016-07-19 Thread Steven Yen
I recently updated my R and RStudio to the latest version and now the 
binary option in the "build" command in devtools stops working.

I went around and used the binary=F option which worked by I get the 
.tar.gz file instead of the .zip file which I prefer.

Does anyone understand the following error message:

status 127
running 'zip' failed

===
setwd("A:/R/yenlib/"); library(devtools)
#build("yenlib",binary=T) # Thisfailed with an error message
build("yenlib",binary=F) # This works

 > build("yenlib",binary=T)
"C:/PROGRA~1/R/R-33~1.1/bin/x64/R" --no-site-file  \
   --no-environ --no-save --no-restore --quiet CMD INSTALL \
   "A:\R\yenlib\yenlib" --build

* installing to library 
'C:/Users/syen01/AppData/Local/Temp/Rtmp8A7KEw/temp_libpath4074149a528e'
* installing *source* package 'yenlib' ...
** R
** data
** preparing package for lazy loading
** help
*** installing help indices
** building package indices
** testing if installed package can be loaded
*** arch - i386
*** arch - x64
* MD5 sums
Warning: running command '"zip" -r9Xq "A:/R/yenlib/yenlib_16.3.zip" 
yenlib' had status 127
running 'zip' failed
* DONE (yenlib)
[1] "A:/R/yenlib/yenlib_16.3.zip"
 >


[[alternative HTML version deleted]]

__
R-help@r-project.org mailing list -- To UNSUBSCRIBE and more, see
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] Label multiple-line plot

2016-07-19 Thread Greg Snow
with 72 lines to label it will be crowded whatever you do.

Here is one option (though I am showing with fewer lines):

x <- sapply(1:15, function(i) cumsum(rnorm(100)))
par(mar=c(5,4,4,3)+0.1)
matplot(1:100, x, type='l', lty=1)
mtext(LETTERS[1:15], side=4, at=x[100,], las=1, line=1)

One way to spread out the labels is:

library(TeachingDemos)
new.y <- spread.labs(x[100,], mindiff = 1.2*strheight('A'),
 min=min(x), max=max(x))
mtext(LETTERS[1:15], side=4, at=new.y, las=1, line=2)


There is also another label spreading function in the plotrix package
and another option for finding space in a plot for labeling lines in
either the rms or Hmisc package.  There are probably other packages
with tools available now as well, but those are the ones that I am
familiar with.



On Tue, Jul 19, 2016 at 1:36 PM, John Wasige  wrote:
> Thanks! It worked  for me.
> matplot(for_jhon$ID, for_jhon[,2:73], type='l')
>
> Any I dea on how I can label multiple-line plot based on column names?
>
> Thanks for your help
>
> John
>
> On Tue, Jul 19, 2016 at 8:46 PM, Greg Snow <538...@gmail.com> wrote:
>>
>> Most attachments get stripped off, so your data did not make it through.
>>
>> But try:
>>
>> matplot(for_jhon$ID, for_jhon[,2:73], type='l')
>>
>>
>> On Tue, Jul 19, 2016 at 12:24 PM, John Wasige 
>> wrote:
>> > Dear all,
>> >
>> > This is to kindly request for your help. I would like to plot my data.
>> >
>> > The R script below gives some plot that is not clear. How can I get a
>> > clear
>> > multiple-line plot. The data is attached herewith.
>> >
>> > ##R Script
>> > for_jhon = read.csv("C:/LVM_share/for_ jhon.csv", header=TRUE, sep=";")
>> > matplot(for_jhon$ID, cbind(for_jhon[,2:73]))
>> >
>> > Thanks for your help
>> >
>> > John
>> > __
>> > R-help@r-project.org mailing list -- To UNSUBSCRIBE and more, see
>> > 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.
>>
>>
>>
>> --
>> Gregory (Greg) L. Snow Ph.D.
>> 538...@gmail.com
>
>
>
>
> --
> John Wasige
> "There are no REGRATES in LIFE, just lessons (Jennifer Aniston)”



-- 
Gregory (Greg) L. Snow Ph.D.
538...@gmail.com

__
R-help@r-project.org mailing list -- To UNSUBSCRIBE and more, see
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] Label multiple-line plot

2016-07-19 Thread John Wasige
Thanks! It worked  for me.
​matplot(for_jhon$ID, for_jhon[,2:73], type='l')

Any I dea on how I can label multiple-line plot based on column names?

Thanks for your help

John

On Tue, Jul 19, 2016 at 8:46 PM, Greg Snow <538...@gmail.com> wrote:

> Most attachments get stripped off, so your data did not make it through.
>
> But try:
>
> ​​
> matplot(for_jhon$ID, for_jhon[,2:73], type='l')
>
>
> On Tue, Jul 19, 2016 at 12:24 PM, John Wasige 
> wrote:
> > Dear all,
> >
> > This is to kindly request for your help. I would like to plot my data.
> >
> > The R script below gives some plot that is not clear. How can I get a
> clear
> > multiple-line plot. The data is attached herewith.
> >
> > ##R Script
> > for_jhon = read.csv("C:/LVM_share/for_ jhon.csv", header=TRUE, sep=";")
> > matplot(for_jhon$ID, cbind(for_jhon[,2:73]))
> >
> > Thanks for your help
> >
> > John
> > __
> > R-help@r-project.org mailing list -- To UNSUBSCRIBE and more, see
> > 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.
>
>
>
> --
> Gregory (Greg) L. Snow Ph.D.
> 538...@gmail.com
>



-- 
John Wasige
"There are no REGRATES in LIFE, just lessons (Jennifer Aniston)”

[[alternative HTML version deleted]]

__
R-help@r-project.org mailing list -- To UNSUBSCRIBE and more, see
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] Sampe numbers

2016-07-19 Thread Frederic Ntirenganya
The problem is that there
are some zeros  while the numbers should range from 1 to  100

As the replace= TRUE if the smallest candidate is 1 may be replaced
because if I can have 1, fifty times this is equal to 50.
Remember that it is a sampling with a replacement

I am still wondering how I can remove the zero value.


Frederic Ntirenganya
Maseno University,
African Maths Initiative,
Kenya.
Mobile:(+254)718492836
Email: fr...@aims.ac.za
https://sites.google.com/a/aims.ac.za/fredo/

On Tue, Jul 19, 2016 at 9:11 PM,  wrote:

> N <- 100
> C <- 50
> x <- numeric(N)
> for (i in 1:N){
>   x[i] <- sample(C-sum(x),1)
> }
> x
> sum(x)
>
>
>
>
> *Frederic Ntirenganya >*
> Sent by: "R-help" 
>
> 07/19/2016 01:41 PM
> To
> "r-help@r-project.org" ,
>
> cc
> Subject
> [R] Sampe numbers
>
>
>
>
> Hi Guys,
> I am trying to sample 100 numbers from 1:100
> i did it  like this:
>
> sample(1:100,100, replace= TRUE)
> but i want again include a constraint that their sum must be equal to 50
>
> How could I do it?
>
> Cheers
>
> Frederic Ntirenganya
> Maseno University,
> African Maths Initiative,
> Kenya.
> Mobile:(+254)718492836
> Email: fr...@aims.ac.za
> https://sites.google.com/a/aims.ac.za/fredo/
>
>
> [[alternative HTML version deleted]]
>
> __
> R-help@r-project.org mailing list -- To UNSUBSCRIBE and more, see
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide
> http://www.R-project.org/posting-guide.html
> 
> and provide commented, minimal, self-contained, reproducible code.
>
>

[[alternative HTML version deleted]]

__
R-help@r-project.org mailing list -- To UNSUBSCRIBE and more, see
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] multiple-line plot

2016-07-19 Thread Greg Snow
Most attachments get stripped off, so your data did not make it through.

But try:

matplot(for_jhon$ID, for_jhon[,2:73], type='l')


On Tue, Jul 19, 2016 at 12:24 PM, John Wasige  wrote:
> Dear all,
>
> This is to kindly request for your help. I would like to plot my data.
>
> The R script below gives some plot that is not clear. How can I get a clear
> multiple-line plot. The data is attached herewith.
>
> ##R Script
> for_jhon = read.csv("C:/LVM_share/for_ jhon.csv", header=TRUE, sep=";")
> matplot(for_jhon$ID, cbind(for_jhon[,2:73]))
>
> Thanks for your help
>
> John
> __
> R-help@r-project.org mailing list -- To UNSUBSCRIBE and more, see
> 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.



-- 
Gregory (Greg) L. Snow Ph.D.
538...@gmail.com

__
R-help@r-project.org mailing list -- To UNSUBSCRIBE and more, see
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] multiple-line plot

2016-07-19 Thread John Wasige
​Dear all,

This is to kindly request for your help. I would like to plot my data.

The R script below gives some plot that is not clear. How can I get a clear
multiple-line plot. The data is attached herewith.

##R Script
for_jhon = read.csv("C:/LVM_share/for_ jhon.csv", header=TRUE, sep=";")
matplot(for_jhon$ID, cbind(for_jhon[,2:73]))​

Thanks for your help

John
__
R-help@r-project.org mailing list -- To UNSUBSCRIBE and more, see
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] pairs: adjusting margins and labeling axes

2016-07-19 Thread Greg Snow
If you want square plots on a rectangular plotting region, then where
do you want the extra space to go?

One option would be to add outer margins to use up the extra space.
The calculations to figure out exactly how much space to put in the
outer margins will probably not be trivial.

Another option would be to not use `pairs`, but use the `layout`
function directly and loops to do your plots (and use the `respect`
argument to `layout`).

On Tue, Jul 19, 2016 at 11:29 AM, michael young
 wrote:
> The default shape for this correlation scatterplot is rectangle.  I changed
> it to square, but then the x-axis spacing between squares are off.  Is
> there an easy way to change x-axis spacing between squares to that of the
> y-axis spacing size?
>
> I decided to hide the name values of the diagonal squares.  I want them
> along the x and y axis instead, outside of the fixed number scale I have.
> I haven't seen any online example of 'pairs' with this and all my searches
> have yielded nothing.  Any ideas?  Thanks
>
> par(pty="s")
> panel.cor <- function(x, y, digits = 2, prefix="", cex.cor, ...)
> {
> usr <- par("usr"); on.exit(par(usr))
> par(usr = c(0, 1, 0, 1),xlog=FALSE,ylog=FALSE)
> # correlation coefficient
> r <- cor(x, y)
> txt <- format(c(r, 0.123456789), digits = digits)[1]
> txt <- paste("r= ", txt, sep = "")
> if(missing(cex.cor)) cex.cor <- 0.8/strwidth(txt)
> text(0.5, 0.6, txt, cex=cex.cor * r)
>
> # p-value calculation
> p <- cor.test(x, y)$p.value
> txt2 <- format(c(p, 0.123456789), digits = digits)[1]
> txt2 <- paste("p= ", txt2, sep = "")
> if(p<0.01) txt2 <- paste("p= ", "<0.01", sep = "")
> text(0.5, 0.4, txt2)
> }
>
> pairs(iris, upper.panel = panel.cor,xlim=c(0.1,10),
> ylim=c(0.1,10),log="xy",text.panel = NULL,pch=".")
>
> [[alternative HTML version deleted]]
>
> __
> R-help@r-project.org mailing list -- To UNSUBSCRIBE and more, see
> 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.



-- 
Gregory (Greg) L. Snow Ph.D.
538...@gmail.com

__
R-help@r-project.org mailing list -- To UNSUBSCRIBE and more, see
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] Sampe numbers

2016-07-19 Thread JLucke
N <- 100
C <- 50
x <- numeric(N)
for (i in 1:N){
  x[i] <- sample(C-sum(x),1)
}
x
sum(x)





Frederic Ntirenganya  
Sent by: "R-help" 
07/19/2016 01:41 PM

To
"r-help@r-project.org" , 
cc

Subject
[R] Sampe numbers






Hi Guys,
I am trying to sample 100 numbers from 1:100
i did it  like this:

sample(1:100,100, replace= TRUE)
but i want again include a constraint that their sum must be equal to 50

How could I do it?

Cheers

Frederic Ntirenganya
Maseno University,
African Maths Initiative,
Kenya.
Mobile:(+254)718492836
Email: fr...@aims.ac.za
https://sites.google.com/a/aims.ac.za/fredo/

 [[alternative HTML version deleted]]

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


[[alternative HTML version deleted]]

__
R-help@r-project.org mailing list -- To UNSUBSCRIBE and more, see
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] Sampe numbers

2016-07-19 Thread Greg Snow
I think that you need to reconsider your conditions.

The smallest number in your candidate set is 1, so if you sample 100
1's they will add to 100 which is greater than 50.  So to have a set
of numbers that sums to 50 you will need to either include negative
numbers, 0's, or sample fewer than 50 values.

If you really meant sum to 500 or sample 10 numbers (still selecting
any single number over 50 will make the sum to large) then this
becomes possible.  There are a few ways to do it, which is best
depends on what you are really trying to accomplish (and note that
with the constraint, the values will not be iid).

One option is rejection sampling:  take a sample, if it sums correctly
then you are done, if not then throw out that sample and try a new
one.  This will be very inefficient, but is easy enough to code and
with a fast enough computer may be acceptable.

Another option is adjustment:  take a sample and compute the sum, if
it is to low then add the difference to one of your values (lowest,
chosen at random, other) or split it up and add to multiple values.
If the sum is too high, then subtract from your numbers until the sum
is correct (checking that you don't go below any lower bounds).

Another option is to round values generated from a Dirichlet
distribution multiplied by your total (you may need to round one value
differently than the default).

Another option is to consider it a balls and urns problem:  You have
500 balls (assuming the true sum is to be 500) that you want
distributed into 100 urns.  If each urn needs at least 1 ball (minimum
value of 1) then put 1 ball in each urn, then for the remaining 400
balls, randomly choose an urn to put each into.  Then count the number
of balls in each urn and those are your 100 numbers that sum to 500
(as long as no balls bounced off the lip of the urn and rolled under
the couch).

There are probably other ways, but these are some things to get you started.

On Tue, Jul 19, 2016 at 11:41 AM, Frederic Ntirenganya
 wrote:
> Hi Guys,
> I am trying to sample 100 numbers from 1:100
> i did it  like this:
>
> sample(1:100,100, replace= TRUE)
> but i want again include a constraint that their sum must be equal to 50
>
> How could I do it?
>
> Cheers
>
> Frederic Ntirenganya
> Maseno University,
> African Maths Initiative,
> Kenya.
> Mobile:(+254)718492836
> Email: fr...@aims.ac.za
> https://sites.google.com/a/aims.ac.za/fredo/
>
> [[alternative HTML version deleted]]
>
> __
> R-help@r-project.org mailing list -- To UNSUBSCRIBE and more, see
> 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.



-- 
Gregory (Greg) L. Snow Ph.D.
538...@gmail.com

__
R-help@r-project.org mailing list -- To UNSUBSCRIBE and more, see
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] pairs: adjusting margins and labeling axes

2016-07-19 Thread michael young
The default shape for this correlation scatterplot is rectangle.  I changed
it to square, but then the x-axis spacing between squares are off.  Is
there an easy way to change x-axis spacing between squares to that of the
y-axis spacing size?

I decided to hide the name values of the diagonal squares.  I want them
along the x and y axis instead, outside of the fixed number scale I have.
I haven't seen any online example of 'pairs' with this and all my searches
have yielded nothing.  Any ideas?  Thanks

par(pty="s")
panel.cor <- function(x, y, digits = 2, prefix="", cex.cor, ...)
{
usr <- par("usr"); on.exit(par(usr))
par(usr = c(0, 1, 0, 1),xlog=FALSE,ylog=FALSE)
# correlation coefficient
r <- cor(x, y)
txt <- format(c(r, 0.123456789), digits = digits)[1]
txt <- paste("r= ", txt, sep = "")
if(missing(cex.cor)) cex.cor <- 0.8/strwidth(txt)
text(0.5, 0.6, txt, cex=cex.cor * r)

# p-value calculation
p <- cor.test(x, y)$p.value
txt2 <- format(c(p, 0.123456789), digits = digits)[1]
txt2 <- paste("p= ", txt2, sep = "")
if(p<0.01) txt2 <- paste("p= ", "<0.01", sep = "")
text(0.5, 0.4, txt2)
}

pairs(iris, upper.panel = panel.cor,xlim=c(0.1,10),
ylim=c(0.1,10),log="xy",text.panel = NULL,pch=".")

[[alternative HTML version deleted]]

__
R-help@r-project.org mailing list -- To UNSUBSCRIBE and more, see
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] Sampe numbers

2016-07-19 Thread Frederic Ntirenganya
Hi Guys,
I am trying to sample 100 numbers from 1:100
i did it  like this:

sample(1:100,100, replace= TRUE)
but i want again include a constraint that their sum must be equal to 50

How could I do it?

Cheers

Frederic Ntirenganya
Maseno University,
African Maths Initiative,
Kenya.
Mobile:(+254)718492836
Email: fr...@aims.ac.za
https://sites.google.com/a/aims.ac.za/fredo/

[[alternative HTML version deleted]]

__
R-help@r-project.org mailing list -- To UNSUBSCRIBE and more, see
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] Concatenate two lists replacing elements with the same name.

2016-07-19 Thread Gabor Grothendieck
Try this:

Reduce(modifyList, list(x, y, z))

On Tue, Jul 19, 2016 at 12:34 PM, Luca Cerone  wrote:
> Dear all,
> I would like to know if there is a function to concatenate two lists
> while replacing elements with the same name.
>
> For example:
>
> x <- list(a=1,b=2,c=3)
> y <- list( b=4, d=5)
> z <- list(a = 6, b = 8, e= 7)
>
> I am looking for a function "concatfun" so that
>
> u <- concatfun(x,y,z)
>
> returns:
>
> u$a=6
> u$b=8
> u$c=3
> u$d=5
> u$e=7
>
> I.e. it combines the 3 lists, but when names have the same value it
> keeps the most recent one.
>
> Does such a function exists?
>
> Thanks for the help,
>
> Cheers,
> Luca
>
> __
> R-help@r-project.org mailing list -- To UNSUBSCRIBE and more, see
> 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.



-- 
Statistics & Software Consulting
GKX Group, GKX Associates Inc.
tel: 1-877-GKX-GROUP
email: ggrothendieck at gmail.com

__
R-help@r-project.org mailing list -- To UNSUBSCRIBE and more, see
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] Concatenate two lists replacing elements with the same name.

2016-07-19 Thread William Dunlap via R-help
concatfun <- function(...) {
   Reduce(f=function(a,b){ a[names(b)] <- b ; a },x=list(...), init=list())
}


Bill Dunlap
TIBCO Software
wdunlap tibco.com

On Tue, Jul 19, 2016 at 9:34 AM, Luca Cerone  wrote:

> Dear all,
> I would like to know if there is a function to concatenate two lists
> while replacing elements with the same name.
>
> For example:
>
> x <- list(a=1,b=2,c=3)
> y <- list( b=4, d=5)
> z <- list(a = 6, b = 8, e= 7)
>
> I am looking for a function "concatfun" so that
>
> u <- concatfun(x,y,z)
>
> returns:
>
> u$a=6
> u$b=8
> u$c=3
> u$d=5
> u$e=7
>
> I.e. it combines the 3 lists, but when names have the same value it
> keeps the most recent one.
>
> Does such a function exists?
>
> Thanks for the help,
>
> Cheers,
> Luca
>
> __
> R-help@r-project.org mailing list -- To UNSUBSCRIBE and more, see
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide
> http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
>

[[alternative HTML version deleted]]

__
R-help@r-project.org mailing list -- To UNSUBSCRIBE and more, see
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] Concatenate two lists replacing elements with the same name.

2016-07-19 Thread Luca Cerone
Dear all,
I would like to know if there is a function to concatenate two lists
while replacing elements with the same name.

For example:

x <- list(a=1,b=2,c=3)
y <- list( b=4, d=5)
z <- list(a = 6, b = 8, e= 7)

I am looking for a function "concatfun" so that

u <- concatfun(x,y,z)

returns:

u$a=6
u$b=8
u$c=3
u$d=5
u$e=7

I.e. it combines the 3 lists, but when names have the same value it
keeps the most recent one.

Does such a function exists?

Thanks for the help,

Cheers,
Luca

__
R-help@r-project.org mailing list -- To UNSUBSCRIBE and more, see
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] Missing rows anova

2016-07-19 Thread Michael Dewey
Presumably it disappears because there is a unique value of ID for eac 
combination of S*x1 so they are indistinguishable.


On 19/07/2016 12:53, Justin Thong wrote:

Why does the S:x1 column disappear (presumably S:x1 goes into ID but I dont
know why)? S is a factor, x1 is a covariate and ID is a factor.

rich.side<-aov(y~S*x1+ID)
summary(rich.side)

Below is the model frame

model.frame(~S*x1+ID)

S x1  ID
1   1 12   A
2   1 12   A
3   1 12   A
4   1 12   A
5   1  0   B
6   1  0   B
7   1  0   B
8   1  0   B
9   1  0   C
10  1  0   C
11  1  0   C
12  1  0   C
13  1  0   D
14  1  0   D
15  1  0   D
16  1  0   D
17  1  0   E
18  1  0   E
19  1  0   E
20  1  0   E
21  1  0   F
22  1  0   F
23  1  0   F
24  1  0   F
25  2  6  AB
26  2  6  AB
27  2  6  AB
28  2  6  AB
29  2  6  AC
30  2  6  AC
31  2  6  AC
32  2  6  AC
33  2  6  AD
34  2  6  AD
35  2  6  AD
36  2  6  AD
37  2  6  AE
38  2  6  AE
39  2  6  AE
40  2  6  AE
41  2  6  AF
42  2  6  AF
43  2  6  AF
44  2  6  AF
45  2  0  BC
46  2  0  BC
47  2  0  BC
48  2  0  BC
49  2  0  BD
50  2  0  BD
51  2  0  BD
52  2  0  BD
53  2  0  BE
54  2  0  BE
55  2  0  BE
56  2  0  BE
57  2  0  BF
58  2  0  BF
59  2  0  BF
60  2  0  BF
61  2  0  CD
62  2  0  CD
63  2  0  CD
64  2  0  CD
65  2  0  CE
66  2  0  CE
67  2  0  CE
68  2  0  CE
69  2  0  CF
70  2  0  CF
71  2  0  CF
72  2  0  CF
73  2  0  DE
74  2  0  DE
75  2  0  DE
76  2  0  DE
77  2  0  DF
78  2  0  DF
79  2  0  DF
80  2  0  DF
81  2  0  EF
82  2  0  EF
83  2  0  EF
84  2  0  EF
85  3  4 ABC
86  3  4 ABC
87  3  4 ABC
88  3  4 ABC
89  3  4 ABD
90  3  4 ABD
91  3  4 ABD
92  3  4 ABD
93  3  4 ABE
94  3  4 ABE
95  3  4 ABE
96  3  4 ABE
97  3  4 ABF
98  3  4 ABF
99  3  4 ABF
100 3  4 ABF
101 3  4 ACD
102 3  4 ACD
103 3  4 ACD
104 3  4 ACD
105 3  4 ACE
106 3  4 ACE
107 3  4 ACE
108 3  4 ACE
109 3  4 ACF
110 3  4 ACF
111 3  4 ACF
112 3  4 ACF
113 3  4 ADE
114 3  4 ADE
115 3  4 ADE
116 3  4 ADE
117 3  4 ADF
118 3  4 ADF
119 3  4 ADF
120 3  4 ADF
121 3  4 AEF
122 3  4 AEF
123 3  4 AEF
124 3  4 AEF
125 3  0 BCD
126 3  0 BCD
127 3  0 BCD
128 3  0 BCD
129 3  0 BCE
130 3  0 BCE
131 3  0 BCE
132 3  0 BCE
133 3  0 BCF
134 3  0 BCF
135 3  0 BCF
136 3  0 BCF
137 3  0 BDE
138 3  0 BDE
139 3  0 BDE
140 3  0 BDE
141 3  0 BDF
142 3  0 BDF
143 3  0 BDF
144 3  0 BDF
145 3  0 BEF
146 3  0 BEF
147 3  0 BEF
148 3  0 BEF
149 3  0 CDE
150 3  0 CDE
151 3  0 CDE
152 3  0 CDE
153 3  0 CDF
154 3  0 CDF
155 3  0 CDF
156 3  0 CDF
157 3  0 CEF
158 3  0 CEF
159 3  0 CEF
160 3  0 CEF
161 3  0 DEF
162 3  0 DEF
163 3  0 DEF
164 3  0 DEF



--
Michael
http://www.dewey.myzen.co.uk/home.html

__
R-help@r-project.org mailing list -- To UNSUBSCRIBE and more, see
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] Regression with ARMA residual

2016-07-19 Thread Mangalani Peter Makananisa

Hi All,

From the library forecast I have fitted a regression model with ARMA residuals 
on a transformed variable diff(log(Y),1). 

What "code(s)" must I use to get the fitted and forecasted values on level 
values ( or original scale of Y) without doing my own manual manipulation?

Please advice

Kind regards

Peter
082 456 4669

-Original Message-
From: Ismail SEZEN [mailto:sezenism...@gmail.com]
Sent: 13 July 2016 04:55 PM
To: Mangalani Peter Makananisa
Cc: r-help@r-project.org
Subject: Re: [R] Dates in R (Year Month)

You can not convert numeric vectors directly to yearmon object. You must 
convert the X variable to character and add “-“ between year and month. Then 
as.yearmon function will work properly.

Please, read help pages of ?as.character, ?strptime and ?as.yearmon.

Example:

> library(zoo)
> aaa <- as.yearmon(c("2015-02”, “2014-06"))
> class(aaa)
[1] "yearmon"

> On 13 Jul 2016, at 16:25, Mangalani Peter Makananisa 
>  wrote:
> 
> Hi All,
> 
> I am trying to convert the vector below to dates please assist I have tried 
> to use information on the links you sent, but it is not working.
> 
> X  = c(201501, 201502, 201503, 201505, 201506, 201507, 201508, 201509, 
> 201510, 201511, 201512, 201601, 201602, 201603, 201604, 201605,
> 201606)
> 
> library(chron, zoo)
> Z = as.yearmon(X)  # it is not working
> 
> please assist
> 
> Kind regards
> Peter
> 
> Please Note: This email and its contents are subject to our email 
> legal notice which can be viewed at 
> http://www.sars.gov.za/Pages/Email-disclaimer.aspx
> 
>   [[alternative HTML version deleted]]
> 
> __
> R-help@r-project.org mailing list -- To UNSUBSCRIBE and more, see 
> 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.

Please Note: This email and its contents are subject to our email legal notice 
which can be viewed at http://www.sars.gov.za/Pages/Email-disclaimer.aspx

[[alternative HTML version deleted]]

__
R-help@r-project.org mailing list -- To UNSUBSCRIBE and more, see
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] frequency of items

2016-07-19 Thread sri vathsan
Hi,

I have a data frame like below.
11,15,12,25
11,12
15,25
134,45,56
46
45,56
15,12
66,45,56,24,14,11,25,12,134

I want to identify the frequency of pairs/triplets or higher that occurs in
the data. Say for example, in above data the occurrence of pairs looks like
below

item No of occurrence
11,12 3
11,25 2
15,12 2
15,25 2
.
.
45,56  3
134,45,562

 and so on

I am trying to write R code for the above and I am finding difficulty to
approach this. Looking forward some help.

Thanks!

-- 

Regards,
Srivathsan.K

[[alternative HTML version deleted]]

__
R-help@r-project.org mailing list -- To UNSUBSCRIBE and more, see
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-es] Pregunta sobre boxplots

2016-07-19 Thread Gerardo Gold
Gracias, Carlos. Lo que quiero es representar de manera clara que los datos
recolectados cada 24 horas no se distribuyen normalmente y que la
dispersion es muy grande, sobre todo al principio del experimento. Tambien
he pensado en hacer un grafico de puntos. La audiencia esperada para esos
graficos no tiene un buen background en estadistica.

La estructura es:

0 hr  Var1   Var2  Var3 . . . . .
24 hrVar1   Var2  Var3 . . . .

Gerardo

2016-07-18 17:44 GMT-05:00 Carlos Ortega :

> Hola,
>
> Cuando dices que es la "mejor representación que he encontrado", ¿a qué te
> refieres?.
> ¿A representar las frecuencias de los datos que consigues?...
>
> Si no estoy entendiendo mal tus datos tienen esta pinta:
>
> Fecha_hora_1 valor_replica_1_1 valor_replica_1_2 ... valor_replica_1_n
> Fecha_hora_2 valor_replica_2_1 valor_replica_2_2 ... valor_replica_2_n
> 
> Fecha_hora_m valor_replica_m_1 valor_replica_m_2 ... valor_replica_m_n
>
> ¿Es así?...
> Si es así, con este tipo de estructura puedes hacer muchos tipos de
> análisis al margen del categórico.
>
> Pero antes de especular, lo mejor es que confirmes la estructura de tus
> datos o envíes un ejemplo para ya sugerirte alguna opción en particular.
>
> Gracias,
> Carlos Ortega
> www.qualityexcellence.es
>
>
> El 18 de julio de 2016, 18:25, Gerardo Gold 
> escribió:
>
>> Colegas:
>>
>> Tengo una pregunta rara, y ofrezco mis disculpas si es una tonteria.
>>
>>  Estoy tratando de hacer unos graficos con datos recolectados cada 24
>> horas, con varias replicas por cada periodo de tiempo. La mejor
>> representacion que he encontrado es hacer boxplots haciendo el tiempo una
>> variable categorica. Sin embargo, me gustaria tener una linea de regresion,
>> o ajuste usando splines o LOESS a travez de las medianas, pero como tengo
>> que hacer el tiempo una variable categorica no puedo hacer el ajuste.
>>
>> Si esto suena demasiado disparatado, podrian por favor sugerir otras
>> alternativas?
>>
>> Gracias anticipadas por su ayuda.
>>
>> Saludos cordiales,
>>
>> Gerardo
>>
>
>
>
> --
> Saludos,
> Carlos Ortega
> www.qualityexcellence.es
>

[[alternative HTML version deleted]]

___
R-help-es mailing list
R-help-es@r-project.org
https://stat.ethz.ch/mailman/listinfo/r-help-es


Re: [R] Has R recently made performance improvements in accumulation?

2016-07-19 Thread boB Rudis
Ideally, you would use a more functional programming approach:

minimal <- function(rows, cols){
  x <- matrix(NA_integer_, ncol = cols, nrow = 0)
  for (i in seq_len(rows)){
x <- rbind(x, rep(i, 10))
  }
  x
}

minimaly <- function(rows, cols){
  x <- matrix(NA_integer_, ncol = cols, nrow = 0)
  do.call(rbind, lapply(seq_len(rows), rep, cols))
}

identical(minimal(100, 100), minimaly(100, 100))
# [1] TRUE

microbenchmark(
  .for=minimal(100, 100),
  .lap=minimaly(100, 100)
)

## Unit: microseconds
##  expr minlq  mean   median   uq  max neval cld
## .for 943.936 1062.3710 1416.1399 1120.259 1366.860 4655.322   100   b
## .lap 111.566  118.1945  160.9058  124.520  146.991 2862.391   100  a

On Tue, Jul 19, 2016 at 10:27 AM, Thierry Onkelinx
 wrote:
> Dear Brent,
>
> I can confirm your timings with
>
> library(microbenchmark)
> microbenchmark(
>   mkFrameForLoop(100, 10),
>   mkFrameForLoop(200, 10),
>   mkFrameForLoop(400, 10)
> )
>
> but profiling your code shown that rbind only uses a small fraction on the
> cpu time used by the function.
>
> profvis::profvis({mkFrameForLoop(100, 10)})
>
> So I cleaned your example further into the function below. Now rbind is
> using most of cpu time. And the timings indicate an O(n^2) relation.
>
> minimal <- function(rows, cols){
>   x <- matrix(NA_integer_, ncol = cols, nrow = 0)
>   for (i in seq_len(rows)){
> x <- rbind(x, rep(i, 10))
>   }
> }
>
> profvis::profvis({minimal(1000, 100)})
>
> timing <- microbenchmark(
>   X50 = minimal(50, 50),
>   X100 = minimal(100, 50),
>   X200 = minimal(200, 50),
>   X400 = minimal(400, 50),
>   X800 = minimal(800, 50),
>   X1600 = minimal(1600, 50)
> )
> timing
> Unit: microseconds
>   exprmin lqmean median  uqmax
> neval cld
>X50199.006212.278233.8444235.728247.3770296.987
>   100 a
>   X100565.693593.957827.8733618.835640.1925   2950.139
>   100 a
>   X200   1804.059   1876.390   2166.1106   1903.370   1938.7115   4263.967
>   100 a
>   X400   6453.913   8755.848   8546.4339   8890.884   8961.7535  13259.024
>   100 a
>   X800  30575.048  32913.186  36555.0118  33093.243  34620.5895 178740.765
>   100  b
>  X1600 130976.429 133674.679 151494.6492 135197.087 137327.1235 292291.385
>   100   c
> timing$N <- as.integer(gsub("X", "", levels(timing$expr)))[timing$expr]
> model <- lm(time ~ poly(N, 4), data = timing)
> summary(model)
>
> Call:
> lm(formula = time ~ poly(N, 4), data = timing)
>
> Residuals:
>   Min1QMedian3Q   Max
> -20518162  -3378940   -130815-45881 142183951
>
> Coefficients:
>   Estimate Std. Error t value Pr(>|t|)
> (Intercept)   33303987 843350  39.490   <2e-16 ***
> poly(N, 4)1 1286962014   20657783  62.299   <2e-16 ***
> poly(N, 4)2  338770077   20657783  16.399   <2e-16 ***
> poly(N, 4)3 222734   20657783   0.0110.991
> poly(N, 4)4   -2260902   20657783  -0.1090.913
> ---
> Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
>
> Residual standard error: 2066 on 595 degrees of freedom
> Multiple R-squared:  0.8746, Adjusted R-squared:  0.8738
> F-statistic:  1038 on 4 and 595 DF,  p-value: < 2.2e-16
> newdata <- data.frame(N = pretty(timing$N, 40))
> newdata$time <- predict(model, newdata = newdata)
> plot(newdata$N, newdata$time, type = "l")
> plot(newdata$N, sqrt(newdata$time), type = "l")
>
> model2 <- lm(sqrt(time) ~ poly(N, 4), data = timing)
> summary(model2)
> Call:
> lm(formula = sqrt(time) ~ poly(N, 4), data = timing)
>
> Residuals:
>Min 1Q Median 3QMax
> -756.3 -202.8  -54.7   -5.5 7416.5
>
> Coefficients:
>  Estimate Std. Error t value Pr(>|t|)
> (Intercept)   3980.36  33.13 120.160  < 2e-16 ***
> poly(N, 4)1 100395.40 811.41 123.730  < 2e-16 ***
> poly(N, 4)2   2191.34 811.41   2.701  0.00712 **
> poly(N, 4)3   -803.54 811.41  -0.990  0.32243
> poly(N, 4)4 82.09 811.41   0.101  0.91945
> ---
> Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
>
> Residual standard error: 811.4 on 595 degrees of freedom
> Multiple R-squared:  0.9626, Adjusted R-squared:  0.9624
> F-statistic:  3829 on 4 and 595 DF,  p-value: < 2.2e-16
>
>
> ir. Thierry Onkelinx
> Instituut voor natuur- en bosonderzoek / Research Institute for Nature and
> Forest
> team Biometrie & Kwaliteitszorg / team Biometrics & Quality Assurance
> Kliniekstraat 25
> 1070 Anderlecht
> Belgium
>
> To call in the statistician after the experiment is done may be no more
> than asking him to perform a post-mortem examination: he may be able to say
> what the experiment died of. ~ Sir Ronald Aylmer Fisher
> The plural of anecdote is not data. ~ Roger Brinner
> The combination of some data and an aching desire for an answer does not
> ensure that a reasonable answer can be extracted from a given body of data.
> ~ John Tukey
>
> 2016-07-19 15:40 GMT+02:00 Brent via 

Re: [R] Has R recently made performance improvements in accumulation?

2016-07-19 Thread Thierry Onkelinx
Dear Brent,

I can confirm your timings with

library(microbenchmark)
microbenchmark(
  mkFrameForLoop(100, 10),
  mkFrameForLoop(200, 10),
  mkFrameForLoop(400, 10)
)

but profiling your code shown that rbind only uses a small fraction on the
cpu time used by the function.

profvis::profvis({mkFrameForLoop(100, 10)})

So I cleaned your example further into the function below. Now rbind is
using most of cpu time. And the timings indicate an O(n^2) relation.

minimal <- function(rows, cols){
  x <- matrix(NA_integer_, ncol = cols, nrow = 0)
  for (i in seq_len(rows)){
x <- rbind(x, rep(i, 10))
  }
}

profvis::profvis({minimal(1000, 100)})

timing <- microbenchmark(
  X50 = minimal(50, 50),
  X100 = minimal(100, 50),
  X200 = minimal(200, 50),
  X400 = minimal(400, 50),
  X800 = minimal(800, 50),
  X1600 = minimal(1600, 50)
)
timing
Unit: microseconds
  exprmin lqmean median  uqmax
neval cld
   X50199.006212.278233.8444235.728247.3770296.987
  100 a
  X100565.693593.957827.8733618.835640.1925   2950.139
  100 a
  X200   1804.059   1876.390   2166.1106   1903.370   1938.7115   4263.967
  100 a
  X400   6453.913   8755.848   8546.4339   8890.884   8961.7535  13259.024
  100 a
  X800  30575.048  32913.186  36555.0118  33093.243  34620.5895 178740.765
  100  b
 X1600 130976.429 133674.679 151494.6492 135197.087 137327.1235 292291.385
  100   c
timing$N <- as.integer(gsub("X", "", levels(timing$expr)))[timing$expr]
model <- lm(time ~ poly(N, 4), data = timing)
summary(model)

Call:
lm(formula = time ~ poly(N, 4), data = timing)

Residuals:
  Min1QMedian3Q   Max
-20518162  -3378940   -130815-45881 142183951

Coefficients:
  Estimate Std. Error t value Pr(>|t|)
(Intercept)   33303987 843350  39.490   <2e-16 ***
poly(N, 4)1 1286962014   20657783  62.299   <2e-16 ***
poly(N, 4)2  338770077   20657783  16.399   <2e-16 ***
poly(N, 4)3 222734   20657783   0.0110.991
poly(N, 4)4   -2260902   20657783  -0.1090.913
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Residual standard error: 2066 on 595 degrees of freedom
Multiple R-squared:  0.8746, Adjusted R-squared:  0.8738
F-statistic:  1038 on 4 and 595 DF,  p-value: < 2.2e-16
newdata <- data.frame(N = pretty(timing$N, 40))
newdata$time <- predict(model, newdata = newdata)
plot(newdata$N, newdata$time, type = "l")
plot(newdata$N, sqrt(newdata$time), type = "l")

model2 <- lm(sqrt(time) ~ poly(N, 4), data = timing)
summary(model2)
Call:
lm(formula = sqrt(time) ~ poly(N, 4), data = timing)

Residuals:
   Min 1Q Median 3QMax
-756.3 -202.8  -54.7   -5.5 7416.5

Coefficients:
 Estimate Std. Error t value Pr(>|t|)
(Intercept)   3980.36  33.13 120.160  < 2e-16 ***
poly(N, 4)1 100395.40 811.41 123.730  < 2e-16 ***
poly(N, 4)2   2191.34 811.41   2.701  0.00712 **
poly(N, 4)3   -803.54 811.41  -0.990  0.32243
poly(N, 4)4 82.09 811.41   0.101  0.91945
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Residual standard error: 811.4 on 595 degrees of freedom
Multiple R-squared:  0.9626, Adjusted R-squared:  0.9624
F-statistic:  3829 on 4 and 595 DF,  p-value: < 2.2e-16


ir. Thierry Onkelinx
Instituut voor natuur- en bosonderzoek / Research Institute for Nature and
Forest
team Biometrie & Kwaliteitszorg / team Biometrics & Quality Assurance
Kliniekstraat 25
1070 Anderlecht
Belgium

To call in the statistician after the experiment is done may be no more
than asking him to perform a post-mortem examination: he may be able to say
what the experiment died of. ~ Sir Ronald Aylmer Fisher
The plural of anecdote is not data. ~ Roger Brinner
The combination of some data and an aching desire for an answer does not
ensure that a reasonable answer can be extracted from a given body of data.
~ John Tukey

2016-07-19 15:40 GMT+02:00 Brent via R-help :

> Subtitle: or, more likely, am I benchmarking wrong?
>
> I am new to R, but I have read that R is a hive of performance pitfalls.
> A very common case is if you try to accumulate results into a typical
> immutable R data structure.
>
> Exhibit A for this confusion is this StackOverflow question on an
> algorithm for O(1) performance in appending to a list:
>
> https://stackoverflow.com/questions/2436688/append-an-object-to-a-list-in-r-in-amortized-constant-time-o1
> The original question was asked in 2010-03-13, and responses have trickled
> in for at least the next 5 years.
>
> Before mentioning my problem, I instead offer an example from someone
> vastly better than me in R, which I am sure should be free of mistakes.
> Consider this interesting article on efficient accumulation in R:
> http://www.win-vector.com/blog/2015/07/efficient-accumulation-in-r/
>
> In that article's first example, it claims that the "obvious “for-loop”
> solution" (accumulation into a data.frame using rbind) 

[R] Has R recently made performance improvements in accumulation?

2016-07-19 Thread Brent via R-help
Subtitle: or, more likely, am I benchmarking wrong?

I am new to R, but I have read that R is a hive of performance pitfalls.  A 
very common case is if you try to accumulate results into a typical immutable R 
data structure.

Exhibit A for this confusion is this StackOverflow question on an algorithm for 
O(1) performance in appending to a list:

https://stackoverflow.com/questions/2436688/append-an-object-to-a-list-in-r-in-amortized-constant-time-o1
The original question was asked in 2010-03-13, and responses have trickled in 
for at least the next 5 years.

Before mentioning my problem, I instead offer an example from someone vastly 
better than me in R, which I am sure should be free of mistakes.  Consider this 
interesting article on efficient accumulation in R:
http://www.win-vector.com/blog/2015/07/efficient-accumulation-in-r/

In that article's first example, it claims that the "obvious “for-loop” 
solution" (accumulation into a data.frame using rbind) is:
is incredibly slow.
...
we pay the cost of processing n*(n+1)/2 rows of data

In other words, what looks like an O(n) algorithm is actually O(n^2).

Sounds awful.  But when I tried executing their example on my machine, I see 
O(n) NOT O(n^2) behavior!

Here is the complete code that I executed:

mkFrameForLoop <- function(nRow,nCol) {
d <- c()
for(i in seq_len(nRow)) {
ri <- mkRow(nCol)
di <- data.frame(ri, stringsAsFactors=FALSE)
d <- rbind(d,di)
}
d
}

mkRow <- function(nCol) {
x <- as.list(rnorm(nCol))
# make row mixed types by changing first column to string
x[[1]] <- ifelse(x[[1]]>0,'pos','neg')
names(x) <- paste('x',seq_len(nCol),sep='.')
x
}

t1 = Sys.time()
x = mkFrameForLoop(100, 10)
tail(x)
t2 = Sys.time()
t2 - t1

t1 = Sys.time()
x = mkFrameForLoop(200, 10)
tail(x)
t2 = Sys.time()
t2 - t1

t1 = Sys.time()
x = mkFrameForLoop(400, 10)
tail(x)
t2 = Sys.time()
t2 - t1

And here is what I got for the execution times:

Time difference of 0.113005876541138 secs
Time difference of 0.205012083053589 secs
Time difference of 0.390022993087769 secs

That's a linear, not quadratic, increase in execution time as a function of 
nRow!  It is NOT what I was expecting, altho it was nice to see.

Notes:
--the functions above are copy and pastes from the article
--the execution time measurement code is all that I added
--yes, that time measurement code is crude
--I am aware of R benchmarking frameworks, in fact, I first tried using 
some of them
--but when I got unexpected results, I went back to the simplest code 
possible, which is the above
--it turns out that I get the same results regardless of how I measure
--each measurement doubles the number of rows (from 100 to 200 to 400), 
which is what should be increasing to bring out rbind's allegedly bad behavior
--my machine has a Xeon E3 1245, 8 GB of RAM, running Win 7 Pro 64 bit, 
using R 3.3.1 (latest release as of today)

Given those unexpected results, you now understand the title of my post.  So, 
has R's runtime somehow gotten greater intelligence recently so that it knows 
when it does not need to copy a data structure?  Maybe in the case above, it 
does a quick shallow copy instead of a deep copy?  Or, am I benchmarking wrong?

What motivated my investigation is that I want to store a bunch of Monte Carlo 
simulation results in a numeric matrix.  When I am finished with my code, I 
know that it will be a massive matrix (e.g. 1,000,000 or more rows).  So I was 
concerned about performance, and went about benchmarking.  Let me know if you 
want to see that benchmark code.  I found that assignment to a matrix does not 
seem to generate a copy, even tho assignment is a mutating operation, so I was 
worried that it might.  But when I investigated deeper, such as with the above 
code, I got worried that I cannot trust my results.


I look forward to any feedback that you can give.

I would especially appreciate any links that explain how you can determine what 
the R runtime actually does.

__
R-help@r-project.org mailing list -- To UNSUBSCRIBE and more, see
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] [FORGED] intersection of two polygons which are not shapefiles

2016-07-19 Thread Adrienne Wootten
Thanks! That's just what I was looking for!
A

On Mon, Jul 18, 2016 at 5:56 PM, Rolf Turner 
wrote:

> On 19/07/16 01:16, Adrienne Wootten wrote:
>
>> All,
>>
>> Greetings! I hope things are going well for all!  I apologize if someone's
>> already answered this and I'm just not finding it.  Here's what I'd like
>> to
>> do, but I'm hitting a brick wall with it.
>>
>> I have two sets of points for which I've already determined which ones
>> points for the boundaries with the chull function.  What I need for what
>> I'm doing is the coordinates where the two resulting polygons overlap.
>> These are not raster grids, nor shapefiles.  What I have to work with are
>> the two data frames with the points at the vertices of each polygon
>> (included below).
>>
>> chx1
>>>
>>   x  y
>> 1  0.5822569 -0.878
>> 2  0.5338428 -0.5883604
>> 3 -0.3442943 -0.6701115
>> 4  -0.7409293  0.2286962
>> 5  0.2147221  0.8061485
>> 6  0.4914146  0.4941865
>> 7  0.5822569 -0.878
>>
>> chx2
>>>
>>   x  y
>> 1  0.7163506 -0.4357497
>> 2  0.6513128 -0.5395180
>> 3   0.1818315 -0.6317423
>> 4  -0.6394281 -0.5765610
>> 5 -0.6044681  0.1831627
>> 6 -0.5799485  0.3231473
>> 7  0.2248476  0.9601908
>> 8  0.7163506 -0.4357497
>>
>>
>> If anyone has any ideas on how to get what I'm after I'd appreciate it!
>> I've tried a lot of things from the raster, rgeos, and more.  Knowing me
>> it's something obvious I'm just not seeing right now.
>>
>
> You can do this very easily using the spatstat package:
>
> library(spatstat)
> X1 <- read.table(textConnection("
>  x  y
> 1  0.5822569 -0.878
> 2  0.5338428 -0.5883604
> 3 -0.3442943 -0.6701115
> 4  -0.7409293  0.2286962
> 5  0.2147221  0.8061485
> 6  0.4914146  0.4941865
> 7  0.5822569 -0.878"))
>
> X2 <- read.table(textConnection("
>x  y
> 1  0.7163506 -0.4357497
> 2  0.6513128 -0.5395180
> 3   0.1818315 -0.6317423
> 4  -0.6394281 -0.5765610
> 5 -0.6044681  0.1831627
> 6 -0.5799485  0.3231473
> 7  0.2248476  0.9601908
> 8  0.7163506 -0.4357497"))
>
> X1 <- reverse.xypolygon(X1) # Vertices are in the wrong
> X2 <- reverse.xypolygon(X2) # (clockwise) order.
> W1 <- owin(poly=X1)
> W2 <- owin(poly=X2)
> WI <- intersect.owin(W1,W2)
>
> plot(union.owin(W1,W2),col="blue",main="")
> plot(WI,add=TRUE,col="red")
>
> HTH
>
> cheers,
>
> Rolf Turner
>
> P.S.  To extract the coordinates of the intersection polygon you can do:
>
> WI$bdry[[1]]
>
> R. T.
>
> --
> Technical Editor ANZJS
> Department of Statistics
> University of Auckland
> Phone: +64-9-373-7599 ext. 88276
>



-- 
Adrienne Wootten
Ph.D Candidate / Graduate Research Assistant
State Climate Office of North Carolina
Department of Marine, Earth and Atmospheric Sciences
North Carolina State University

[[alternative HTML version deleted]]

__
R-help@r-project.org mailing list -- To UNSUBSCRIBE and more, see
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] Anova() type iii SS plots and diagnostics

2016-07-19 Thread Fox, John
Dear Pamela,

I'm afraid that this question seems confused. Type III sums of squares are 
computed for a linear model -- it's the *linear model* that you'd check, in the 
usual manner, for outliers, etc., and this has nothing to do with how the ANOVA 
table for the model is computed.

Best,
 John

-
John Fox, Professor
McMaster University
Hamilton, Ontario
Canada L8S 4M4
Web: socserv.mcmaster.ca/jfox



> -Original Message-
> From: R-help [mailto:r-help-boun...@r-project.org] On Behalf Of Pamela
> Wong
> Sent: July 18, 2016 10:39 PM
> To: r-help@r-project.org
> Subject: [R] Anova() type iii SS plots and diagnostics
> 
> I am wondering if there is a way to plot results and model diagnostics (to 
> check
> for outliers, homoscedasticity, normality, collinearity) using type III sums 
> of
> squares in R __
> R-help@r-project.org mailing list -- To UNSUBSCRIBE and more, see
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide http://www.R-project.org/posting-
> guide.html
> and provide commented, minimal, self-contained, reproducible code.

__
R-help@r-project.org mailing list -- To UNSUBSCRIBE and more, see
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] Summing up layers

2016-07-19 Thread maettuw
Hello everyone

I have gridded output data from a climate model. I would like to take the 
average of these grids and plot them. I have tried to plot it with gglot or 
image.plot (field/raster package). I have to do this extensively meaning that a 
manual way is not what I am looking for. My data has three dimensions. 
Longitude,Latitude and the variable itself. Below is a code that works 
perfectly fine. This is a plot for one single year.

H3 <- ((detrended_TOA_rad[,,28:28]))
image.plot(lon,lat,H3,col=tim.colors())

setwd("/scratch/maettu/maettu/data/fup013/echam5/net_Radiation")
spatial_TOA_net_Rad <- open.nc("spatial_yearmean_TS_TOA_net_radiation.nc")
TOA_net_rad <- read.nc(spatial_TOA_net_Rad)
detrended_TOA_rad <- TOA_net_rad$ttr
lon_TOA <- TOA_net_rad$lon
lat_TOA <- TOA_net_rad$lat


Down here I have tried to make a plot over several years, which did not work!!

H3 <- ((detrended_TOA_rad[,,28:38]))
image.plot(lon,lat,H3,col=tim.colors())

-->

Error in seq.default(zlim[1], zlim[2], , nlevel) :
  'from' cannot be NA, NaN or infinite

I googled that error message but didnt get smart out of it.
The ggplot way worked better but the problem is that I cannot create a dataset 
that contains the coordinates as well as the values. Therefore the projection 
is not correct.

melted1 <- melt(H3)

ggplot(aes(x = Var1, y = Var2, fill = value), data = melted1) +geom_tile() 
+geom_tile() + scale_fill_distiller(palette = "Set1") +
  geom_polygon(data=map,aes(x=long, y=lat, group=group), colour="black", 
fill="transparent")

I have also tried to aggregate with the "overlay" function from the raster 
packege. However I cannot coerce my array files into a raster...
So the main problem is that I cannot aggrregate the layers.

Appreciate the help! Matthias

[[alternative HTML version deleted]]

__
R-help@r-project.org mailing list -- To UNSUBSCRIBE and more, see
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] Missing rows anova

2016-07-19 Thread Justin Thong
Why does the S:x1 column disappear (presumably S:x1 goes into ID but I dont
know why)? S is a factor, x1 is a covariate and ID is a factor.

rich.side<-aov(y~S*x1+ID)
summary(rich.side)

Below is the model frame

model.frame(~S*x1+ID)

S x1  ID
1   1 12   A
2   1 12   A
3   1 12   A
4   1 12   A
5   1  0   B
6   1  0   B
7   1  0   B
8   1  0   B
9   1  0   C
10  1  0   C
11  1  0   C
12  1  0   C
13  1  0   D
14  1  0   D
15  1  0   D
16  1  0   D
17  1  0   E
18  1  0   E
19  1  0   E
20  1  0   E
21  1  0   F
22  1  0   F
23  1  0   F
24  1  0   F
25  2  6  AB
26  2  6  AB
27  2  6  AB
28  2  6  AB
29  2  6  AC
30  2  6  AC
31  2  6  AC
32  2  6  AC
33  2  6  AD
34  2  6  AD
35  2  6  AD
36  2  6  AD
37  2  6  AE
38  2  6  AE
39  2  6  AE
40  2  6  AE
41  2  6  AF
42  2  6  AF
43  2  6  AF
44  2  6  AF
45  2  0  BC
46  2  0  BC
47  2  0  BC
48  2  0  BC
49  2  0  BD
50  2  0  BD
51  2  0  BD
52  2  0  BD
53  2  0  BE
54  2  0  BE
55  2  0  BE
56  2  0  BE
57  2  0  BF
58  2  0  BF
59  2  0  BF
60  2  0  BF
61  2  0  CD
62  2  0  CD
63  2  0  CD
64  2  0  CD
65  2  0  CE
66  2  0  CE
67  2  0  CE
68  2  0  CE
69  2  0  CF
70  2  0  CF
71  2  0  CF
72  2  0  CF
73  2  0  DE
74  2  0  DE
75  2  0  DE
76  2  0  DE
77  2  0  DF
78  2  0  DF
79  2  0  DF
80  2  0  DF
81  2  0  EF
82  2  0  EF
83  2  0  EF
84  2  0  EF
85  3  4 ABC
86  3  4 ABC
87  3  4 ABC
88  3  4 ABC
89  3  4 ABD
90  3  4 ABD
91  3  4 ABD
92  3  4 ABD
93  3  4 ABE
94  3  4 ABE
95  3  4 ABE
96  3  4 ABE
97  3  4 ABF
98  3  4 ABF
99  3  4 ABF
100 3  4 ABF
101 3  4 ACD
102 3  4 ACD
103 3  4 ACD
104 3  4 ACD
105 3  4 ACE
106 3  4 ACE
107 3  4 ACE
108 3  4 ACE
109 3  4 ACF
110 3  4 ACF
111 3  4 ACF
112 3  4 ACF
113 3  4 ADE
114 3  4 ADE
115 3  4 ADE
116 3  4 ADE
117 3  4 ADF
118 3  4 ADF
119 3  4 ADF
120 3  4 ADF
121 3  4 AEF
122 3  4 AEF
123 3  4 AEF
124 3  4 AEF
125 3  0 BCD
126 3  0 BCD
127 3  0 BCD
128 3  0 BCD
129 3  0 BCE
130 3  0 BCE
131 3  0 BCE
132 3  0 BCE
133 3  0 BCF
134 3  0 BCF
135 3  0 BCF
136 3  0 BCF
137 3  0 BDE
138 3  0 BDE
139 3  0 BDE
140 3  0 BDE
141 3  0 BDF
142 3  0 BDF
143 3  0 BDF
144 3  0 BDF
145 3  0 BEF
146 3  0 BEF
147 3  0 BEF
148 3  0 BEF
149 3  0 CDE
150 3  0 CDE
151 3  0 CDE
152 3  0 CDE
153 3  0 CDF
154 3  0 CDF
155 3  0 CDF
156 3  0 CDF
157 3  0 CEF
158 3  0 CEF
159 3  0 CEF
160 3  0 CEF
161 3  0 DEF
162 3  0 DEF
163 3  0 DEF
164 3  0 DEF

-- 
Yours sincerely,
Justin

*I check my email at 9AM and 4PM everyday*
*If you have an EMERGENCY, contact me at +447938674419(UK) or
+60125056192(Malaysia)*

[[alternative HTML version deleted]]

__
R-help@r-project.org mailing list -- To UNSUBSCRIBE and more, see
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] Factors in nlme: singularity in backsolve error

2016-07-19 Thread Alexander Shenkin

Hello all,

I tried posting on r-sig-ME, but it's not going through.  So I thought I 
would try here.  Apologies if this ends up being a cross-post if it 
eventually goes through there


I am parameterizing exponential fits for some metabolic scaling models. 
I have done this in lmer already, without problem, using logged 
dependent and independent variables.  However, I would now like to 
incorporate other parameters that aren't necessarily exponentially 
related to the dependent variable.  Hence, I've turned to nlme, but I 
don't have much experience with it.  Apologies in advance for newbie 
mistakes.


With the code below, I am getting the following error.  I'm guessing 
that it has something to do with the 'site' factor being misspecified:


Error in nlme.formula(dep ~ scaling_fun(alpha, beta, ind, site), data = 
scale_df,  :

  Singularity in backsolve at level 0, block 1

When I fit a simpler function that does not involve 'site', the model 
seems to work correctly.


Any thoughts would be greatly appreciated!  I've added the dput object 
at the bottom of the email to avoid distraction...


Thanks,
Allie


> head(scale_df)
  depind spp site
2  0.28069471 -0.0322841 157A
3 -0.69719050 -1.2568901 183A
4  0.29252012  0.1592420 246A
5  0.72030740 -0.3282789 154A
6 -0.08601891  0.3623756 110A
7  0.30793594  0.2230840 154A


scaling_fun <- function(alpha, beta, ind, site) {
return(beta + ind^alpha + site*(ind^alpha))
}

# results in singularity in backsolve error

nlme(dep ~ trait_scaling_fun(alpha, beta, ind, site),
 data = scale_df,
 fixed = list(alpha + beta + site ~ 1), random = alpha ~ 1|spp,
 start = list(fixed = c(0.7, 0, 1)))


##
# simpler function converges #
##

scaling_fun <- function(alpha, beta, ind) {
return(beta + ind^alpha)
}

nlme(dep ~ scaling_fun(alpha, beta, ind),
 data = scale_df,
 fixed = list(alpha + beta ~ 1), random = alpha ~ 1|spp,
 start = list(fixed = c(0.7, 0)))


# dput for data

scale_df = structure(list(dep = structure(c(0.280694709349609, 
-0.697190504779316,

0.29252012378561, 0.720307403285465, -0.0860189064012942, 0.307935942585907,
0.862281389721553, -0.265859454130292, -0.662753116661401, 
-0.281590648203109,

1.17791119955669, 0.631729125553481, -0.556750637034717, -0.43297610538,
0.78495323332528, 0.355622302250686, -1.26351337875593, 0.343598303486811,
-1.17350365646135, 0.660405635460007, -0.0444581036277384, 
0.722773887452235,
-0.0522058979055465, 0.65594868736652, -0.471798559774723, 
0.580581816047275,
0.514408011465886, -0.736138172672309, 0.70774857676606, 
-0.0810467855114514,

2.87816317587089, -0.22475032263, -0.497119013269447, 0.031618943500335,
1.93466683745019, 1.0320275291612, 0.874248439252774, -0.327339718896358,
-0.854016320390115, -1.14806351731618, -0.380316794335449, 
0.298626272822636,
-0.600084009068698, 0.836832381353764, -1.02141473458357, 
-0.0626961141428343,
-1.28013982373819, -0.126310567347687, -0.998914167194147, 
-0.966356030272785,

-0.0311625866367497, -0.456212263484603, 2.22345644904786, 1.06808914728274,
2.04676618838053, 2.57057502230498, -1.1354565673953, -0.617861267279932,
-1.29433182429494, 0.0191424880072437, -0.31539529665173, 
-0.161941249027648,

-0.582100497448459, 1.62041284412996, 4.12468548155128, 0.181382553309611,
0.395970527486673, 0.949165316445473, 0.536508784640814, 1.88628468601066,
0.460108411105217, 1.50237201805536, -0.248693954731227, 0.604338152411222,
-0.511271387322619, -0.273625967472597, 2.54906943158448, 1.25797640582933,
-0.615499666866161, -0.83022966646613, -0.569798846726899, 1.41607035018281,
0.806876783115155, 3.03500743372039, 0.371311069097266, 0.605523223251687,
0.960528362908707, 0.536545286280021, -1.07783545300095, 0.145069857506269,
-0.390394315762374, -0.524314770428875, -0.205267560118326, 2.8197556390129,
0.548519073742583, 2.25991992849428, -0.274997856331803, -0.748248845217783,
1.27303941400619, 0.23884797601163, -0.582228991340883, -0.380185805876727,
-0.0858712197617979, -0.181441660791883, -0.122600427784157,
1.44111920784939, 0.163470634515893, -0.322754741867578, -1.08355296391174,
-0.27943230825155, -0.367498411981987, -1.25029544932326, 1.96800896885018,
0.741382953698137, 0.525448033325593, 0.58486584900894, -1.39293857038504,
0.0401527810688047, -0.428292900882259, -0.384773565802791, 
-0.347537789997483,

0.259461537066455, 1.04953062690743, 0.905652281045827, 0.921052829745329,
0.52450347921673, 1.71595979136317, -0.474652211889487, 0.986966207710898,
-0.256654313404628, 2.99291703386829, -0.874752003654369, 1.01178420482112,
5.77381325416311, -0.167933435302975, -0.746365585369451, 2.80467857407177,
1.95866740218257, -0.503053147492697, -0.0160574809046825, 
0.0121332489691898,

0.155574669334708, 0.916506565829058, -0.448074058381429, 0.316834873904141,
1.04981848557227, -0.824609748750358, 1.06744876285501, 0.063863664699137,

Re: [R] ... distribution based on mean values and quantiles in R?

2016-07-19 Thread Martin Maechler
> Jim Lemon 
> on Tue, 19 Jul 2016 16:20:49 +1000 writes:

> Hi Daniel,
> Judging by the numbers you mention, the distribution is either very
> skewed or not at all normal. If you look at this:

> plot(c(0,0.012,0.015,0.057,0.07),c(0,0.05,0.4,0.05,0),type="b")

> you will see the general shape of whatever distribution produced these
> summary statistics. Did the paper give any hints as to what the model
> distribution might be?

Yes, that's the correct question:  At first, it's not about
plotting, about *fitting* a distribution with the desired
properties, and then you can easily visualize it.

So, what were the data?  If they are 'amounts' of any kind, they
are necessarily non-negative often always positive, and hence
--- according to John W Tukey --- should be analyzed after
taking logs (Tukey's "first aid transformation" for *any*
amounts data).

Taking logs, and analyzing means to consider a normal
("Gaussian") distribution for the log()  which is
equivalent to fitting a  lognormal distribution -- R functions [dpqrr]lnorm() 
to the original data.  I'd strongly recommend doing that.

And I did so, finding out, however that if indeed it is the
*mean* and the 15% and 95% quantiles,  a log normal is not
fitting.  Here is the reproducible R code .. with lots of
comments :


##
## MM  strongly recommends to fit a  log-normal distribution .. however it does 
*not* fit

## The "data statistics"
qlN <- c(q15 = 0.012, q85 = 0.057) # Quantiles original scale
mlN <- 0.015

(qn <- log(qlN))  # Quantiles log scale [assumed normal 
(Gaussian) here]
## as the Gaussian is symmetri, the mid value of the two quantiles is the mean 
and median :
(mu <- mean(qn))   # -3.644
(medlN <- exp(mu)) # 0.02615
## an estimate for the median(.)  -- but  it is *larger* than the mean = 0.015 
above !
## ===> Log-normal does *NOT* fit well :

## From the help page, we learn that
##E[lN] = exp(mu + 1/2 sigma^2){and it is trivial that}
##   median[lN] = exp(mu)
## where here, our medLn is a (moment) estimate of median[lN]

## If the number were different, this would solve the problem :
## Consequently, a (moment / plugin) estimate for sigma is
(e12sig <- mlN / medlN) ## ~= exp( 1/2 sigma^2)
(sig2 <- 2*log(e12sig)) ## ~=  sigma^2  [--- is *NEGATIVE* (<==> 'est.median' > 
mean !)]
(sig <- sqrt(sig2)) ## ~=  sigma[here of course 'NaN' with a warning !]


My conclusion would be that other distributions (than the
log-normal; the normal  is out of question !!) have to be
considered, if you want to be serious about the topic.

Maybe the poweRlaw package (https://cloud.r-project.org/package=poweRlaw)
may help you (it has 4 vignettes, the last being a nice JSS publication).

The above is a "cute" non-standard problem in any case: to fit very skewed
distributions, given two quantiles and the mean only, and the
approach taken by the "poweRlawyers", namely to minimize the KS
(Kolmogorov-Smirnoff) decrepancy seems a good start to me.

Martin Maechler,
ETH Zurich



> Jim


> On Tue, Jul 19, 2016 at 7:11 AM, gcchenleidenuniv
>  wrote:
>> Hi all,
>> 
>> I need to draw density curves based on some published data. But in the 
article only mean value (0.015 ) and quantiles (Q0.15=0.012 , Q0.85=0.057) were 
given. So I was thinking if it is possible to plot density curves solely based 
on the mean value and quantiles. The dnorm(x, mean, sd, log) function needs the 
standard deviation which was not mentioned, so it can not be used in this 
situation.
>> 
>> Many thanks!!
>> Daniel
>> [[alternative HTML version deleted]]
>> 
>> __
>> R-help@r-project.org mailing list -- To UNSUBSCRIBE and more, see
>> https://stat.ethz.ch/mailman/listinfo/r-help
>> PLEASE do read the posting guide 
http://www.R-project.org/posting-guide.html
>> and provide commented, minimal, self-contained, reproducible code.

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

__
R-help@r-project.org mailing list -- To UNSUBSCRIBE and more, see
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] update of formula with conditional rhs introduces unwantedparentheses

2016-07-19 Thread Gerrit Eichner

Dear list members,

is anyone aware of a (simple) strategy to prevent update() from
introducing parentheses around a conditional rhs of a formula
when updating, e.g., only its lhs? Simple example:

update( y ~ a|b, log(.) ~ .)

gives log(y) ~ (a | b), but I want/need log(y) ~ a | b


I do know how to extract various components of a formula and how to
"re-build" formulae from such single components, but I would like to
avoid that.

 Thx for any hint  --  Gerrit

-
Dr. Gerrit Eichner   Mathematical Institute, Room 212
gerrit.eich...@math.uni-giessen.de   Justus-Liebig-University Giessen
Tel: +49-(0)641-99-32104  Arndtstr. 2, 35392 Giessen, Germany
Fax: +49-(0)641-99-32109http://www.uni-giessen.de/eichner

__
R-help@r-project.org mailing list -- To UNSUBSCRIBE and more, see
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] R2WinBUGS with Multivariate Logistic Regression

2016-07-19 Thread Uwe Ligges



On 17.07.2016 05:33, Christopher Kelvin via R-help wrote:

Dear R-User,
I have written a simple code to analyze some data using Bayesian logistic 
regression via the R2WinBUGS package. The code when run in WinBUGS stops 
WinBUGS from running it and using the package returns no results also.



I'd suggest to reduce it to a WinBUGS only problem, try OpenBUGS (which 
iw more recent) and if it still fails, ask on the WinBUGS/OpenBUGS 
mailing list. A trap in WinBUGS is not an R related problem.


Best,
Uwe Ligges





I attach herewith, the code and a sample of the dataset.

Any suggestion will be greatly appreciated.

Chris Guure
Biostatistics Department

University of Ghana




library(R2WinBUGS)
library(coda)model1<-function(){
for (i in 1:N) {
# likelihood function
ms[i] ~ dbin( p [ i ], N )
logit(p [ i ] ) <- alpha + bage*age[ i ] + bpam*pam[i ] + bpah*pah[i]

}
### prior for intercept
alpha ~ dnorm(0,0.0001)

# prior for slopes
bage ~ dnorm(0,0.0001)
bpam ~ dnorm(0,0.0001)
bpah ~ dnorm(0,0.0001)


# OR for alpha
or.age<-exp(bage)
# OR for hbp
or.pam <- exp(bpam)
# OR for fdm
or.pah <- exp(bpah)

}

data=cbind(ms=c(0, 0, 0, 0, 1, 0, 0, 0, 0, 0, 0, 0, 0, 1, 0, 0, 0, 0, 0, 0),
age=c(77, 83, 75, 78, 75, 83, 85, 80, 80, 85, 76, 77, 80, 76, 88, 77, 81, 78, 
85, 81),
pam=c(0, 0, 1, 1, 0, 1, 0, 0, 0, 0, 1, 0, 0, 0, 0, 0, 0, 0, 0, 1),
pah=c(1, 1, 0, 0, 1, 0, 1, 1, 1, 1, 0, 0, 1, 1, 1, 1, 1, 1, 0, 0))
N=20

 initial values 
###
lineinits <- function(){
list(alpha=1,bage= 0.05, bpam=0.05, bpah=0.05)
}

## CODA output
lineout1 <- bugs(data, lineinits, c("bage", "alpha","bpam", "bpah", "or.age", "or.pam", 
"or.pah"), model1,n.iter = 11000, n.burnin = 1000, n.chains = 2, codaPkg = T, DIC = TRUE)


### Posterior summaries
line.coda <- read.bugs(lineout1)
summary(line.coda)

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



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


Re: [R] How to plot density distribution based on mean values and quantiles in R?

2016-07-19 Thread Jim Lemon
Hi Daniel,
Judging by the numbers you mention, the distribution is either very
skewed or not at all normal. If you look at this:

plot(c(0,0.012,0.015,0.057,0.07),c(0,0.05,0.4,0.05,0),type="b")

you will see the general shape of whatever distribution produced these
summary statistics. Did the paper give any hints as to what the model
distribution might be?

Jim


On Tue, Jul 19, 2016 at 7:11 AM, gcchenleidenuniv
 wrote:
> Hi all,
>
> I need to draw density curves based on some published data. But in the 
> article only mean value (0.015 ) and quantiles (Q0.15=0.012 , Q0.85=0.057) 
> were given. So I was thinking if it is possible to plot density curves solely 
> based on the mean value and quantiles. The dnorm(x, mean, sd, log) function 
> needs the standard deviation which was not mentioned, so it can not be used 
> in this situation.
>
> Many thanks!!
> Daniel
> [[alternative HTML version deleted]]
>
> __
> R-help@r-project.org mailing list -- To UNSUBSCRIBE and more, see
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.

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


[R] How to plot density distribution based on mean values and quantiles in R?

2016-07-19 Thread gcchenleidenuniv
Hi all,

I need to draw density curves based on some published data. But in the article 
only mean value (0.015 ) and quantiles (Q0.15=0.012 , Q0.85=0.057) were given. 
So I was thinking if it is possible to plot density curves solely based on the 
mean value and quantiles. The dnorm(x, mean, sd, log) function needs the 
standard deviation which was not mentioned, so it can not be used in this 
situation.

Many thanks!!
Daniel
[[alternative HTML version deleted]]

__
R-help@r-project.org mailing list -- To UNSUBSCRIBE and more, see
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] Anova() type iii SS plots and diagnostics

2016-07-19 Thread Pamela Wong
I am wondering if there is a way to plot results and model diagnostics (to 
check for outliers, homoscedasticity, normality, collinearity) using type III 
sums of squares in R
__
R-help@r-project.org mailing list -- To UNSUBSCRIBE and more, see
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] barplot colour problem

2016-07-19 Thread Abdoulaye Sarr
Thank you Marc,

The typo was causing the problem, solved now.

Regards,
Fipou

On Mon, Jul 18, 2016 at 8:38 PM, Marc Schwartz  wrote:

>
> > On Jul 18, 2016, at 1:06 PM, Abdoulaye Sarr 
> wrote:
> >
> > I am doing a basic bar plot which works but the color of bars positive
> > (green) and negative (brown) don’t show up from the below command:
> >
> > barplot(z, ylim=c(-2,2), col=ifelse(x>0,"brown","green »))
> >
> > any help? or other methods?
> >
> > fipou
>
>
> Presuming that the above is a direct copy and paste, your ifelse()
> statement is using 'x' to determine the color, rather than 'z'. Presumably
> a typo?
>
> This works, for example, with 'z' as a vector:
>
>   z <- seq(from = -5, to = 5)
>
>   barplot(z, col = ifelse(z > 0, "brown", "green"))
>
>
> Regards,
>
> Marc Schwartz
>
>
>

[[alternative HTML version deleted]]

__
R-help@r-project.org mailing list -- To UNSUBSCRIBE and more, see
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.