Re: [R] What exactly is an dgCMatrix-class. There are so many attributes.

[David Winsemius]

> On Oct 20, 2017, at 11:11 AM, C W  wrote:
> 
> Dear R list,
> 
> I came across dgCMatrix. I believe this class is associated with sparse
> matrix.

Yes. See:

 help('dgCMatrix-class', pack=Matrix)

If Martin Maechler happens to respond to this you should listen to him rather 
than anything I write. Much of what the Matrix package does appears to be 
magical to one such as I.

> 
> I see there are 8 attributes to train$data, I am confused why are there so
> many, some are vectors, what do they do?
> 
> Here's the R code:
> 
> library(xgboost)
> data(agaricus.train, package='xgboost')
> data(agaricus.test, package='xgboost')
> train <- agaricus.train
> test <- agaricus.test
> attributes(train$data)
> 

I got a bit of an annoying surprise when I did something similar. It appearred 
to me that I did not need to load the xgboost library since all that was being 
asked was "where is the data" in an object that should be loaded from that 
library using the `data` function. The last command asking for the attributes 
filled up my console with a 100K length vector (actually 2 of such vectors). 
The `str` function returns a more useful result.

> data(agaricus.train, package='xgboost')
> train <- agaricus.train
> names( attributes(train$data) )
[1] "i""p""Dim"  "Dimnames" "x""factors"  "class"   
> str(train$data)
Formal class 'dgCMatrix' [package "Matrix"] with 6 slots
  ..@ i   : int [1:143286] 2 6 8 11 18 20 21 24 28 32 ...
  ..@ p   : int [1:127] 0 369 372 3306 5845 6489 6513 8380 8384 10991 ...
  ..@ Dim : int [1:2] 6513 126
  ..@ Dimnames:List of 2
  .. ..$ : NULL
  .. ..$ : chr [1:126] "cap-shape=bell" "cap-shape=conical" "cap-shape=convex" 
"cap-shape=flat" ...
  ..@ x   : num [1:143286] 1 1 1 1 1 1 1 1 1 1 ...
  ..@ factors : list()

> Where is the data, is it in $p, $i, or $x?

So the "data" (meaning the values of the sparse matrix) are in the @x leaf. The 
values all appear to be the number 1. The @i leaf is the sequence of row 
locations for the values entries while the @p items are somehow connected with 
the columns (I think, since 127 and 126=number of columns from the @Dim leaf 
are only off by 1). 

Doing this > colSums(as.matrix(train$data))
  cap-shape=bellcap-shape=conical 
 3693 
cap-shape=convex   cap-shape=flat 
2934 2539 
   cap-shape=knobbed cap-shape=sunken 
 644   24 
 cap-surface=fibrous  cap-surface=grooves 
18674 
   cap-surface=scaly   cap-surface=smooth 
2607 2035 
 cap-color=brown   cap-color=buff 
1816  
# now snipping the rest of that output.



Now this makes me think that the @p vector gives you the cumulative sum of 
number of items per column:

> all( cumsum( colSums(as.matrix(train$data)) ) == train$data@p[-1] )
[1] TRUE

> 
> Thank you very much!
> 
>   [[alternative HTML version deleted]]

Please read the Posting Guide. Your code was not mangled in this instance, but 
HTML code often arrives in an unreadable mess.

> 
> ______
> 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.

David Winsemius
Alameda, CA, USA

'Any technology distinguishable from magic is insufficiently advanced.'   
-Gehm's Corollary to Clarke's Third Law

__
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] What exactly is an dgCMatrix-class. There are so many attributes.

[David Winsemius]

> On Oct 21, 2017, at 7:50 AM, Martin Maechler  
> wrote:
> 
>>>>>> C W 
>>>>>>on Fri, 20 Oct 2017 15:51:16 -0400 writes:
> 
>> Thank you for your responses.  I guess I don't feel
>> alone. I don't find the documentation go into any detail.
> 
>> I also find it surprising that,
> 
>>> object.size(train$data)
>> 1730904 bytes
> 
>>> object.size(as.matrix(train$data))
>> 6575016 bytes
> 
>> the dgCMatrix actually takes less memory, though it
>> *looks* like the opposite.
> 
> to whom?
> 
> The whole idea of these sparse matrix classes in the 'Matrix'
> package (and everywhere else in applied math, CS, ...) is that
> 1. they need  much less memory   and
> 2. matrix arithmetic with them can be much faster because it is based on
>   sophisticated sparse matrix linear algebra, notably the
>   sparse Cholesky decomposition for solve() etc.
> 
> Of course the efficency only applies if most of the
> matrix entries _are_ 0.
> You can measure the  "sparsity" or rather the  "density", of a
> matrix by
> 
>  nnzero(A) / length(A)
> 
> where length(A) == nrow(A) * ncol(A)  as for regular matrices
> (but it does *not* integer overflow)
> and nnzero(.) is a simple utility from Matrix
> which -- very efficiently for sparseMatrix objects -- gives the
> number of nonzero entries of the matrix.
> 
> All of these classes are formally defined classes and have
> therefore help pages.  Here  ?dgCMatrix-class  which then points
> to  ?CsparseMatrix-class  (and I forget if Rstudio really helps
> you find these ..; in emacs ESS they are found nicely via the usual key)
> 
> To get started, you may further look at  ?Matrix _and_  ?sparseMatrix
> (and possibly the Matrix package vignettes --- though they need
> work -- I'm happy for collaborators there !)
> 
> Bill Dunlap's comment applies indeed:
> In principle all these matrices should work like regular numeric
> matrices, just faster with less memory foot print if they are
> really sparse (and not just formally of a sparseMatrix class)
>  ((and there are quite a few more niceties in the package))
> 
> Martin Maechler
> (here, maintainer of 'Matrix')
> 
> 
>> On Fri, Oct 20, 2017 at 3:22 PM, David Winsemius 
>> wrote:
> 
>>>> On Oct 20, 2017, at 11:11 AM, C W  wrote:
>>>> 
>>>> Dear R list,
>>>> 
>>>> I came across dgCMatrix. I believe this class is associated with sparse
>>>> matrix.
>>> 
>>> Yes. See:
>>> 
>>> help('dgCMatrix-class', pack=Matrix)
>>> 
>>> If Martin Maechler happens to respond to this you should listen to him
>>> rather than anything I write. Much of what the Matrix package does appears
>>> to be magical to one such as I.
>>> 
>>>> 
>>>> I see there are 8 attributes to train$data, I am confused why are there
>>> so
>>>> many, some are vectors, what do they do?
>>>> 
>>>> Here's the R code:
>>>> 
>>>> library(xgboost)
>>>> data(agaricus.train, package='xgboost')
>>>> data(agaricus.test, package='xgboost')
>>>> train <- agaricus.train
>>>> test <- agaricus.test
>>>> attributes(train$data)
>>>> 
>>> 
>>> I got a bit of an annoying surprise when I did something similar. It
>>> appearred to me that I did not need to load the xgboost library since all
>>> that was being asked was "where is the data" in an object that should be
>>> loaded from that library using the `data` function. The last command asking
>>> for the attributes filled up my console with a 100K length vector (actually
>>> 2 of such vectors). The `str` function returns a more useful result.
>>> 
>>>> data(agaricus.train, package='xgboost')
>>>> train <- agaricus.train
>>>> names( attributes(train$data) )
>>> [1] "i""p""Dim"  "Dimnames" "x""factors"
>>> "class"
>>>> str(train$data)
>>> Formal class 'dgCMatrix' [package "Matrix"] with 6 slots
>>> ..@ i   : int [1:143286] 2 6 8 11 18 20 21 24 28 32 ...
>>> ..@ p   : int [1:127] 0 369 372 3306 5845 6489 6513 8380 8384 10991
>>> ...
>>> ..@ Dim : int [1:2] 6513 126
>>> ..@ Dimnames:List of 2
>>> .. ..$ : NULL
>>> .. ..$ : chr [1

Re: [R] Syntax for fit.contrast

[David Winsemius]

> On Oct 22, 2017, at 6:04 AM, Sorkin, John  wrote:
> 
> I have a model (run with glm) that has a factor, type. Type has two levels, 
> "general" and "regional". I am trying to get estimates (and SEs) for the 
> model with type="general" and type ="regional" using fit.contrast

 ?fit.contrast
No documentation for ‘fit.contrast’ in specified packages and libraries:
you could try ‘??fit.contrast’

Perhaps the gmodels function of that name?

> but I can't get the syntax of the coefficients to use in fit.contrast 
> correct. I hope someone can show me how to use fit.contrast, or some other 
> method to get estimate with SEs. (I know I can use the variance co-variance 
> matrix, but I would rather not have to code the linear contrast my self from 
> the coefficients of the matrix)
> 

I'm having trouble understanding what you are trying to extract. There are only 
2 levels so there is really only one interesting contrast ("general" vs 
"regional") , and it's magnitude would be reported by just typing `model`, and 
it's SE would show up in output of `summary(model)`.

I'm thinking you should pick one of the examples in gmodels::fit.contrast that 
most resembles your real problem,  post it,  and  and then explain what 
difficulties you are having with interpretation.

-- 
David.


> Thank  you,
> 
> John
> 
> 
> My model:
> 
> model=glm(events~type,family=poisson(link=log),offset=log(SS),data=data)
> 
> 
> Model details:
> 
>> summary(data$type)
> 
> general regional
>  16   16
> 
>> levels(data$type)
> [1] "general"  "regional"
> 
>> contrasts(data$type)
> regional
> general 0
> regional1
> 
> 
> I have tried the following syntax for fit.contrast
> 
> fit.contrast(model,type,c(1,0))
> and get an error:
> Error in `[[<-`(`*tmp*`, varname, value = cmat) :
>  no such index at level 1
> 
> 
>> fit.contrast(model,type,c(0,1),showall=TRUE)
> and get an error:
> Error in `[[<-`(`*tmp*`, varname, value = cmat) :
>  no such index at level 1
> 
> 
> 
>> fit.contrast(model,type,c(1,-1),showall=TRUE)
> and get an error:
> Error in `[[<-`(`*tmp*`, varname, value = cmat) :
>  no such index at level 1
> 
> 
>> fit.contrast(model,type,c(0))
> and get an error:
> Error in make.contrasts(coeff, ncol(coeff)) :
>  Too many contrasts specified. Must be less than the number of factor levels 
> (columns).
> 
>> fit.contrast(model,type,c(1))
> Error in make.contrasts(coeff, ncol(coeff)) :
> and get an error
>  Too many contrasts specified. Must be less than the number of factor levels 
> (columns).
> 
> 
> 
> 
> 
> 
> 
> 
> John David Sorkin M.D., Ph.D.
> Professor of Medicine
> Chief, Biostatistics and Informatics
> University of Maryland School of Medicine Division of Gerontology and 
> Geriatric Medicine
> Baltimore VA Medical Center
> 10 North Greene Street
> GRECC (BT/18/GR)
> Baltimore, MD 21201-1524
> (Phone) 410-605-7119
> (Fax) 410-605-7913 (Please call phone number above prior to faxing)
> 
> 
>   [[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.

David Winsemius
Alameda, CA, USA

'Any technology distinguishable from magic is insufficiently advanced.'   
-Gehm's Corollary to Clarke's Third Law

__
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] Syntax for fit.contrast (from package gmodels)

[David Winsemius]

> On Oct 22, 2017, at 3:56 PM, Sorkin, John  wrote:
> 
> David,
> Thank you for responding to my post.
> 
> Please consider the following output (typeregional is a factor having two 
> levels, "regional" vs. "general"):
> Call:
> glm(formula = events ~ type, family = poisson(link = log), data = data, 
> offset = log(SS))
> 
> Deviance Residuals: 
> Min   1Q   Median   3Q  Max  
> -43.606  -17.295   -4.6514.204   38.421  
> 
> Coefficients:
>  Estimate Std. Error z value Pr(>|z|)
> (Intercept)  -2.528300.01085 -233.13   <2e-16 ***
> typeregional  0.337880.01641   20.59   <2e-16 ***
> 
> Let's forget for a moment that the model is a Poisson regression and pretend 
> that the output is from a simple linear regression, e.g. from lm.
> 
> To get the estimate for "general" one simply needs to use the value of the 
> intercept i.e. -2.5830. Similarly to get the 95% CI of general one simply 
> needs to compute -2.52830-(1.96*0.01085) and -2.52830+(1.96*0.01085).
> 
> To get the estimate for "regional" one needs to compute intercept + 
> typeregional, i.e. -2.52830 + 0.33788. To get the 95% CI is somewhat more 
> difficult as one needs to use results from the variance-covariance matix, 
> specifically the variance of intercept, the variance of "regional", and the 
> covariance of (intercept,"regional") which involves:
> var =  var(intercept) + var(regional) +2*(covar(intercept,regional)),
> and then get the SE of the variance
> SE=sqrt(var)
> 95% CI = intercept + regional - 1.95*SE and intercept + regional + 1.95*SE.
> 
> I was hoping that a contrast statement could be written that would give me 
> the point estimate and SE for "general" and its SE and another contrast 
> statement could be written that would give me the point estimate and SE for 
> "general" and it SE without my having to work directly with the 
> variance-covariance matrix. I tried doing this using the fit.contrast 
> statements (from the gmodels package):

I'm guessing that the second contrast you were hoping for was actually for 
"regional".

Contrasts, hence the name, are for differences between two levels (or more 
accurately between the means on the scale specified by the link parameter. In 
the absence of another level the only other reference point would be a value of 
zero or perhaps the value you specified by your offset term.

-- 
David


> 
> fit.contrast(model,type,c(1,0),showall=TRUE)
> fit.contrast(model,type,c(0,1),showall=TRUE)
> 
> and received the error message, 
> Error in `[[<-`(`*tmp*`, varname, value = c(0, 1)) : 
>   no such index at level 1
> 
> Perhaps fit.contrast is not the way to accomplish my goal. Perhaps my goal 
> can be accomplished without a contrast statement, but I don't know how.
> 
> Thank you,
> John
>  
> 
> 
> John David Sorkin M.D., Ph.D.
> Professor of Medicine
> Chief, Biostatistics and Informatics
> University of Maryland School of Medicine Division of Gerontology and 
> Geriatric Medicine
> Baltimore VA Medical Center
> 10 North Greene Street
> GRECC (BT/18/GR)
> Baltimore, MD 21201-1524
> (Phone) 410-605-7119
> (Fax) 410-605-7913 (Please call phone number above prior to faxing) 
> 
> 
> 
> From: David Winsemius 
> Sent: Sunday, October 22, 2017 1:20 PM
> To: Sorkin, John
> Cc: r-help@r-project.org
> Subject: Re: [R] Syntax for fit.contrast
>  
> 
> > On Oct 22, 2017, at 6:04 AM, Sorkin, John  wrote:
> > 
> > I have a model (run with glm) that has a factor, type. Type has two levels, 
> > "general" and "regional". I am trying to get estimates (and SEs) for the 
> > model with type="general" and type ="regional" using fit.contrast
> 
>  ?fit.contrast
> No documentation for ‘fit.contrast’ in specified packages and libraries:
> you could try ‘??fit.contrast’
> 
> Perhaps the gmodels function of that name?
> 
> > but I can't get the syntax of the coefficients to use in fit.contrast 
> > correct. I hope someone can show me how to use fit.contrast, or some other 
> > method to get estimate with SEs. (I know I can use the variance co-variance 
> > matrix, but I would rather not have to code the linear contrast my self 
> > from the coefficients of the matrix)
> > 
> 
> I'm having trouble understanding what you are trying to extract. There are 
> only 2 levels so there is really only one interesting contrast ("general" vs 
> "regional") , and it's magnitude would be reported by just typing `model`, 
> and it's SE would show up in

Re: [R] Syntax for fit.contrast (from package gmodels)

[David Winsemius]

> On Oct 22, 2017, at 5:01 PM, Sorkin, John  wrote:
> 
> David,
> Again  you have my thanks!.
> You are correct. What I want is not technically a contrast. What I want is 
> the estimate for "regional" and its SE.

There needs to be a reference value for the contrast. Contrasts are 
differences. I gave you the choice of two references (treatment constrast or 
the offset value you specified). Pick one or suggest an alternate value. 
Typical alternate values are the global mean or zero.

Read ?predict.glm

"se.fit logical switch indicating if standard errors are required."


> I don't mind if I get these on the log scale; I can get the anti-log. Can you 
> suggest  how I can get the point estimate and its SE for "regional"? The 
> predict function will give the point estimate, but not (to my knowledge) the 
> SE.


> Thank you,
> John
> 
> John David Sorkin M.D., Ph.D.
> Professor of Medicine
> Chief, Biostatistics and Informatics
> University of Maryland School of Medicine Division of Gerontology and 
> Geriatric Medicine
> Baltimore VA Medical Center
> 10 North Greene Street
> GRECC (BT/18/GR)
> Baltimore, MD 21201-1524
> (Phone) 410-605-7119
> (Fax) 410-605-7913 (Please call phone number above prior to faxing) 
> 
> 
> 
> From: David Winsemius 
> Sent: Sunday, October 22, 2017 7:56 PM
> To: Sorkin, John
> Cc: r-help@r-project.org
> Subject: Re: [R] Syntax for fit.contrast (from package gmodels)
>  
> 
> > On Oct 22, 2017, at 3:56 PM, Sorkin, John  wrote:
> > 
> > David,
> > Thank you for responding to my post.
> > 
> > Please consider the following output (typeregional is a factor having two 
> > levels, "regional" vs. "general"):
> > Call:
> > glm(formula = events ~ type, family = poisson(link = log), data = data, 
> > offset = log(SS))
> > 
> > Deviance Residuals: 
> > Min   1Q   Median   3Q  Max  
> > -43.606  -17.295   -4.6514.204   38.421  
> > 
> > Coefficients:
> >  Estimate Std. Error z value Pr(>|z|)
> > (Intercept)  -2.528300.01085 -233.13   <2e-16 ***
> > typeregional  0.337880.01641   20.59   <2e-16 ***
> > 
> > Let's forget for a moment that the model is a Poisson regression and 
> > pretend that the output is from a simple linear regression, e.g. from lm.
> > 
> > To get the estimate for "general" one simply needs to use the value of the 
> > intercept i.e. -2.5830. Similarly to get the 95% CI of general one simply 
> > needs to compute -2.52830-(1.96*0.01085) and -2.52830+(1.96*0.01085).
> > 
> > To get the estimate for "regional" one needs to compute intercept + 
> > typeregional, i.e. -2.52830 + 0.33788. To get the 95% CI is somewhat more 
> > difficult as one needs to use results from the variance-covariance matix, 
> > specifically the variance of intercept, the variance of "regional", and the 
> > covariance of (intercept,"regional") which involves:
> > var =  var(intercept) + var(regional) +2*(covar(intercept,regional)),
> > and then get the SE of the variance
> > SE=sqrt(var)
> > 95% CI = intercept + regional - 1.95*SE and intercept + regional + 1.95*SE.
> > 
> > I was hoping that a contrast statement could be written that would give me 
> > the point estimate and SE for "general" and its SE and another contrast 
> > statement could be written that would give me the point estimate and SE for 
> > "general" and it SE without my having to work directly with the 
> > variance-covariance matrix. I tried doing this using the fit.contrast 
> > statements (from the gmodels package):
> 
> I'm guessing that the second contrast you were hoping for was actually for 
> "regional".
> 
> Contrasts, hence the name, are for differences between two levels (or more 
> accurately between the means on the scale specified by the link parameter. In 
> the absence of another level the only other reference point would be a value 
> of zero or perhaps the value you specified by your offset term.
> 
> -- 
> David
> 
> 
> > 
> > fit.contrast(model,type,c(1,0),showall=TRUE)
> > fit.contrast(model,type,c(0,1),showall=TRUE)
> > 
> > and received the error message, 
> > Error in `[[<-`(`*tmp*`, varname, value = c(0, 1)) : 
> >   no such index at level 1
> > 
> > Perhaps fit.contrast is not the way to accomplish my goal. Perhaps my goal 
> > can be accomplished without a contrast statement, but I don't know how.
> > 
> > Thank you,
> > John

Re: [R] Syntax for fit.contrast (from package gmodels)

[David Winsemius]

> On Oct 22, 2017, at 5:26 PM, Sorkin, John  wrote:
> 
> David,
> predict.glm and se.fit were exactly what I was looking for.

The default 'se' delivered for a listed contrast is not for that particular 
level per se, but rather for the difference between that level and the 
non-listed factor level (its mean) which is "embedded" as it were in the 
Intercept. The se for the Intercept is for the mean of all the reference levels 
jointly being zero.

Apologies for any mis-steps in this effort. I'm enjoying my fourth IPA. I 
usually do not drink and derive.

-- 
David.


> Many  thanks!
> John
> 
> John David Sorkin M.D., Ph.D.
> Professor of Medicine
> Chief, Biostatistics and Informatics
> University of Maryland School of Medicine Division of Gerontology and 
> Geriatric Medicine
> Baltimore VA Medical Center
> 10 North Greene Street
> GRECC (BT/18/GR)
> Baltimore, MD 21201-1524
> (Phone) 410-605-7119
> (Fax) 410-605-7913 (Please call phone number above prior to faxing) 
> 
> 
> 
> From: David Winsemius 
> Sent: Sunday, October 22, 2017 8:15 PM
> To: Sorkin, John
> Cc: r-help@r-project.org
> Subject: Re: [R] Syntax for fit.contrast (from package gmodels)
>  
> 
> > On Oct 22, 2017, at 5:01 PM, Sorkin, John  wrote:
> > 
> > David,
> > Again  you have my thanks!.
> > You are correct. What I want is not technically a contrast. What I want is 
> > the estimate for "regional" and its SE.
> 
> There needs to be a reference value for the contrast. Contrasts are 
> differences. I gave you the choice of two references (treatment constrast or 
> the offset value you specified). Pick one or suggest an alternate value. 
> Typical alternate values are the global mean or zero.
> 
> Read ?predict.glm
> 
> "se.fit logical switch indicating if standard errors are required."
> 
> 
> > I don't mind if I get these on the log scale; I can get the anti-log. Can 
> > you suggest  how I can get the point estimate and its SE for "regional"? 
> > The predict function will give the point estimate, but not (to my 
> > knowledge) the SE.
> 
> 
> > Thank you,
> > John
> > 
> > John David Sorkin M.D., Ph.D.
> > Professor of Medicine
> > Chief, Biostatistics and Informatics
> > University of Maryland School of Medicine Division of Gerontology and 
> > Geriatric Medicine
> > Baltimore VA Medical Center
> > 10 North Greene Street
> > GRECC (BT/18/GR)
> > Baltimore, MD 21201-1524
> > (Phone) 410-605-7119
> > (Fax) 410-605-7913 (Please call phone number above prior to faxing) 
> > 
> > 
> > 
> > From: David Winsemius 
> > Sent: Sunday, October 22, 2017 7:56 PM
> > To: Sorkin, John
> > Cc: r-help@r-project.org
> > Subject: Re: [R] Syntax for fit.contrast (from package gmodels)
> >  
> > 
> > > On Oct 22, 2017, at 3:56 PM, Sorkin, John  
> > > wrote:
> > > 
> > > David,
> > > Thank you for responding to my post.
> > > 
> > > Please consider the following output (typeregional is a factor having two 
> > > levels, "regional" vs. "general"):
> > > Call:
> > > glm(formula = events ~ type, family = poisson(link = log), data = data, 
> > > offset = log(SS))
> > > 
> > > Deviance Residuals: 
> > > Min   1Q   Median   3Q  Max  
> > > -43.606  -17.295   -4.6514.204   38.421  
> > > 
> > > Coefficients:
> > >  Estimate Std. Error z value Pr(>|z|)
> > > (Intercept)  -2.528300.01085 -233.13   <2e-16 ***
> > > typeregional  0.337880.01641   20.59   <2e-16 ***
> > > 
> > > Let's forget for a moment that the model is a Poisson regression and 
> > > pretend that the output is from a simple linear regression, e.g. from lm.
> > > 
> > > To get the estimate for "general" one simply needs to use the value of 
> > > the intercept i.e. -2.5830. Similarly to get the 95% CI of general one 
> > > simply needs to compute -2.52830-(1.96*0.01085) and 
> > > -2.52830+(1.96*0.01085).
> > > 
> > > To get the estimate for "regional" one needs to compute intercept + 
> > > typeregional, i.e. -2.52830 + 0.33788. To get the 95% CI is somewhat more 
> > > difficult as one needs to use results from the variance-covariance matix, 
> > > specifically the variance of intercept, the variance of "regional", and 
> > > the covariance of (intercept,"regional") which involves:
&

Re: [R] Problem when trying to run Java in R:

[David Winsemius]

otonmail.com), Swiss-based encrypted email.
>   [[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.

David Winsemius
Alameda, CA, USA

'Any technology distinguishable from magic is insufficiently advanced.'   
-Gehm's Corollary to Clarke's Third Law

__
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] Problem Subsetting Rows that Have NA's

[David Winsemius]

> On Oct 25, 2017, at 6:57 AM, BooBoo  wrote:
> 
> On 10/25/2017 4:38 AM, Ista Zahn wrote:
>> On Tue, Oct 24, 2017 at 3:05 PM, BooBoo  wrote:
>>> This has every appearance of being a bug. If it is not a bug, can someone
>>> tell me what I am asking for when I ask for "x[x[,2]==0,]". Thanks.
>> You are asking for elements of x where the second column is equal to zero.
>> 
>> help("==")
>> 
>> and
>> 
>> help("[")
>> 
>> explain what happens when missing values are involved. I agree that
>> the behavior is surprising, but your first instinct when you discover
>> something surprising should be to read the documentation, not to post
>> to this list. After having read the documentation you may post back
>> here if anything remains unclear.
>> 
>> Best,
>> Ista
>> 
>>>> #here is the toy dataset
>>>> x <- rbind(c(1,1),c(2,2),c(3,3),c(4,0),c(5,0),c(6,NA),
>>> +   c(7,NA),c(8,NA),c(9,NA),c(10,NA)
>>> + )
>>>> x
>>>   [,1] [,2]
>>>  [1,]11
>>>  [2,]22
>>>  [3,]33
>>>  [4,]40
>>>  [5,]50
>>>  [6,]6   NA
>>>  [7,]7   NA
>>>  [8,]8   NA
>>>  [9,]9   NA
>>> [10,]   10   NA
>>>> #it contains rows that have NA's
>>>> x[is.na(x[,2]),]
>>>  [,1] [,2]
>>> [1,]6   NA
>>> [2,]7   NA
>>> [3,]8   NA
>>> [4,]9   NA
>>> [5,]   10   NA
>>>> #seems like an unreasonable answer to a reasonable question
>>>> x[x[,2]==0,]
>>>  [,1] [,2]
>>> [1,]40
>>> [2,]50
>>> [3,]   NA   NA
>>> [4,]   NA   NA
>>> [5,]   NA   NA
>>> [6,]   NA   NA
>>> [7,]   NA   NA
>>>> #this is more what I was expecting
>>>> x[which(x[,2]==0),]
>>>  [,1] [,2]
>>> [1,]40
>>> [2,]50
>>> __
>>> 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.
> 
> I wanted to know if this was a bug so that I could report it if so. You say 
> it is not, so you answered my question. As far as me not reading the 
> documentation, I challenge anyone to read the cited help pages and predict 
> the observed behavior based on the information given in those pages.

Some of us do share (or at least remember feeling) your pain. The ?Extract page 
is long and complex and there are several features that I find non-intuitive. 
But they are deemed desirable by others. I think I needed to read that page 
about ten times (with multiple different problems that needed explication) 
before it started to sink in. You are apparently on that same side of the split 
opinions on the feature of returning rows with logical NA's as I am. I've 
learned to use `which`, and I push back when the conoscienti says it's not 
needed.

 After you read it a few more times you may come to a different opinion. Many 
people come to R with preconceived notions of what words like "equals" or 
"list" or "vector" mean and then complain about the documentation. You would be 
better advised to spend more time studying the language. The help pages are 
precise but terse, and you need to spend time with the examples and with other 
tutorial material to recognize the gotcha's.

Here's a couple of possibly helpful rules regarding "[[" and "[" and logical 
indexing:

Nothing _equals_ NA.
Selection operations with NA logical index item return NA.  (Justified as a 
warning feature as I understand it.)
"[" always returns a list.
"[[" returns only one thing, but even that thing could be a list.
Generally you want "[[" if you plan on testing for equality with a vector.

The "R Inferno" by Burns is an effort to detail many more of the unexpected or 
irregular aspects of R (mostly inherited from S).

-- 
Best of luck in your studies.


> 
> __
> 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.

David Winsemius
Alameda, CA, USA

'Any technology distinguishable from magic is insufficiently advanced.'   
-Gehm's Corollary to Clarke's Third Law

__
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] Problem Subsetting Rows that Have NA's

[David Winsemius]

> On Oct 25, 2017, at 11:17 AM, David Winsemius  wrote:
> 
> 
>> On Oct 25, 2017, at 6:57 AM, BooBoo  wrote:
>> 
>> On 10/25/2017 4:38 AM, Ista Zahn wrote:
>>> On Tue, Oct 24, 2017 at 3:05 PM, BooBoo  wrote:
>>>> This has every appearance of being a bug. If it is not a bug, can someone
>>>> tell me what I am asking for when I ask for "x[x[,2]==0,]". Thanks.
>>> You are asking for elements of x where the second column is equal to zero.
>>> 
>>> help("==")
>>> 
>>> and
>>> 
>>> help("[")
>>> 
>>> explain what happens when missing values are involved. I agree that
>>> the behavior is surprising, but your first instinct when you discover
>>> something surprising should be to read the documentation, not to post
>>> to this list. After having read the documentation you may post back
>>> here if anything remains unclear.
>>> 
>>> Best,
>>> Ista
>>> 
>>>>> #here is the toy dataset
>>>>> x <- rbind(c(1,1),c(2,2),c(3,3),c(4,0),c(5,0),c(6,NA),
>>>> +   c(7,NA),c(8,NA),c(9,NA),c(10,NA)
>>>> + )
>>>>> x
>>>>  [,1] [,2]
>>>> [1,]11
>>>> [2,]22
>>>> [3,]33
>>>> [4,]40
>>>> [5,]50
>>>> [6,]6   NA
>>>> [7,]7   NA
>>>> [8,]8   NA
>>>> [9,]9   NA
>>>> [10,]   10   NA
>>>>> #it contains rows that have NA's
>>>>> x[is.na(x[,2]),]
>>>> [,1] [,2]
>>>> [1,]6   NA
>>>> [2,]7   NA
>>>> [3,]8   NA
>>>> [4,]9   NA
>>>> [5,]   10   NA
>>>>> #seems like an unreasonable answer to a reasonable question
>>>>> x[x[,2]==0,]
>>>> [,1] [,2]
>>>> [1,]40
>>>> [2,]50
>>>> [3,]   NA   NA
>>>> [4,]   NA   NA
>>>> [5,]   NA   NA
>>>> [6,]   NA   NA
>>>> [7,]   NA   NA
>>>>> #this is more what I was expecting
>>>>> x[which(x[,2]==0),]
>>>> [,1] [,2]
>>>> [1,]40
>>>> [2,]50
>>>> __
>>>> 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.
>> 
>> I wanted to know if this was a bug so that I could report it if so. You say 
>> it is not, so you answered my question. As far as me not reading the 
>> documentation, I challenge anyone to read the cited help pages and predict 
>> the observed behavior based on the information given in those pages.
> 
> Some of us do share (or at least remember feeling) your pain. The ?Extract 
> page is long and complex and there are several features that I find 
> non-intuitive. But they are deemed desirable by others. I think I needed to 
> read that page about ten times (with multiple different problems that needed 
> explication) before it started to sink in. You are apparently on that same 
> side of the split opinions on the feature of returning rows with logical NA's 
> as I am. I've learned to use `which`, and I push back when the conoscienti 
> says it's not needed.


horrible misspelling of cognoscenti 


> After you read it a few more times you may come to a different opinion. Many 
> people come to R with preconceived notions of what words like "equals" or 
> "list" or "vector" mean and then complain about the documentation. You would 
> be better advised to spend more time studying the language. The help pages 
> are precise but terse, and you need to spend time with the examples and with 
> other tutorial material to recognize the gotcha's.
> 
> Here's a couple of possibly helpful rules regarding "[[" and "[" and logical 
> indexing:
> 
> Nothing _equals_ NA.
> Selection operations with NA logical index item return NA.  (Justified as a 
> warning feature as I understand it.)
> "[" always returns a list.

That's not true or even half true. "[" always returns a list if it's first 
argument is a list and it only has two arguments.

If X is a list and you ask for X[vector] you get a list

If you ask for X[vector, ] you may get a l

Re: [R] nomogram function error

[David Winsemius]

> On Oct 26, 2017, at 9:42 AM, Annalisa VanderWyden via R-help 
>  wrote:
> 
> 
> Hi R-help,
> 
>  
> 
> I have fit a cox ph model to my data, but have beenreceiving an error when 
> trying to fit a model to the nomogram. Here is the codeand corresponding 
> error:
> 
>  
> 
>  
> 
>> nomogramCF = nomogram(cph_age6_40avp4, 

The only package that I know of having a function named `nomogram` is rms. 
> 
> +lp.at= seq(-10,10,by =0.5),lp = TRUE,
> 
> +  
> 
> +  funlabel="5year survival",
> 
> +  fun=surv5.CFdata5,

That should be a function name being passed. But you haven't shown us any 
function definiton.

> fun.at= c(0.01,seq(0.1,0.9,by=0.2),0.99))
> 
> Error inapprox(fu[s], xseq[s], fat, ties = mean) : 
> 
>   needat least two non-NA values to interpolate
> 
>  
> 
>  
> 
> I have fit similar nomograms based on cox ph models usingsimilar code, so I’m 
> not sure what I’m doing wrong this time.

We are at least as unsure as you and most probably more so. You have not 
offered a reproducible problem.

> I heard there might be some an issue because of a softwareupdate?

There's often a NEWS for packages. Type:

help(pac=rms)

> 
>  
> 
> Thanks,
> 
> Annie
> 
> 
>   [[alternative HTML version deleted]]

Please read the Posting Guide linked to in every rhelp posting. 


-- 
David.
> 
> __
> 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.

David Winsemius
Alameda, CA, USA

'Any technology distinguishable from magic is insufficiently advanced.'   
-Gehm's Corollary to Clarke's Third Law

__
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] 3D Plot with Date Axis?

[David Winsemius]

> On Oct 27, 2017, at 12:22 PM, Alex Restrepo  wrote:
> 
> Hello,
> 
> I would like to format the X axis of the plot created via the scatterplot3d 
> function or any other function which will work.
> 
> Here is an example of what I am trying to do:
> 
> library("scatterplot3d")
> 
> mydf=data.frame(rate=c(1,1,4,4), age=c(2,2,5,5), market_date=c('2017-01-01', 
> '2017-02-02', '2017-03-03', '2017-04-04'))
> 
> scatterplot3d(mydf$market_date, mydf$rate, mydf$age, xaxt="n")
> 
> axis.Date(1, mydf$market_date, format="%Y-%m-%d")
> 
> I tried to hide the x axis, but it looks like xaxt is not working for the 
> scatterplot3d function.

Please read the "?scatterplot3d page more carefully. I think you should be 
looking for the x.ticklabs parameter.

-- 
David.
> 
> Using the above example, could someone please provide me with some guidance 
> and help?
> 
> Many Thanks,
> 
> Alex
> 
> 
>   [[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] Scatterplot3d :: Rotating x tick labels by x degrees

[David Winsemius]

> On Oct 31, 2017, at 8:55 AM, Olivier Crouzet  
> wrote:
> 
> Hi Alex,
> 
> this should be related to the "las" argument of "par()" but
> actually it does not seem to be parametered in scatterplot3d.
> Searching the net for "scatterplot3d las" provides a link to:
> 
> https://stackoverflow.com/questions/25458652/specifying-the-orientation-of-the-axes-labels-in-scatterplot3d
> 
> You may try the solution that is provided in this link or consider using
> alternate packages (like rgl or the plotly packages which one may be
> more powerfull as far as I can judge). However I can't help more. It
> seems ggplot does not produce 3d plots (but it looks like it can
> interact with plotly when using 3d plots).

Olivier;

The cited SO solutions (mine being the first one)  were addressing the request 
for rotated "axis"-labels rather than rotated "axtick"-labels, although the 
general strategy of looking at the code and using the `text`-function with 
xpd=TRUE (rather than the `mtext`-function where needed in the definitions of 
mytext and mytext2) should apply. 

Alex;

I'd encourage you to demonstrate more initiative rather than expecting us to be 
on-call code servants.  I've decided to limit my gratis coding time to 15 
minutes daily. I think this might take me an hour or more. I'm still available 
for on-list code review on this effort. 

As a start, I'd suggest downloading the scatterplot3d package code and then 
open up scatterplot3d.R. Find the comment

## label tick marks

... and perhaps decide whether you need to use `text` rather than `mtext`. 
(`text` can rotate by any amount, 
but may need "xpd" to be set TRUE. `mtext` is limited to 90 degree increments.) 
I've used up my time, so the next move is yours.
-- 
David.


> 
> Olivier.
> 
> On Mon, 30 Oct 2017
> 23:56:02 + Alex Restrepo  wrote:
> 
>> Hi,
>> 
>> I would like to rotate the x axis tick labels by 45 degrees.   Using
>> the code below, could someone please provide an example?   Many
>> Thanks In Advance, Alex
>> 
>> library("scatterplot3d")
>> mydf=data.frame(rate=seq(158, 314)
>>,age=seq(1, 157)
>>,market_date=seq(as.Date("2000/1/1"), as.Date
>> ("2003/1/1"), by="7 days"))
>> 
>> mydf$market_date=as.Date(mydf$market_date, format="%Y-%m-%d")
>> 
>> scatterplot3d(mydf$market_date
>>  ,mydf$rate
>>  ,mydf$age
>>  ,x.ticklabs = seq(as.Date("2000/1/1"), as.Date
>> ("2003/1/1"), by="330 days"))
>> 
>> 
>>  [[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.
> 
> 
> -- 
>  Olivier Crouzet, PhD
>  /Assistant Professor/
>  @LLING - Laboratoire de Linguistique de Nantes
>UMR6310 CNRS / Université de Nantes
>  /Guest Researcher/
>  @UMCG (University Medical Center Groningen)
>ENT department
>Rijksuniversiteit Groningen
> 
> __
> 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.

David Winsemius
Alameda, CA, USA

'Any technology distinguishable from magic is insufficiently advanced.'   
-Gehm's Corollary to Clarke's Third Law

__
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] beta binomial distribution installation

[David Winsemius]

 
> Warning in install.packages :
> 
>  dependency Biobase is not available
> 
> trying URL 
> 'https://cran.rstudio.com/bin/windows/contrib/3.4/TailRank_3.1.3.zip'
> 
> Content type 'application/zip' length 331270 bytes (323 KB)
> 
> downloaded 323 KB
> 
> 
> 
> package TailRank successfully unpacked and MD5 sums checked
> 
> 
> 
> The downloaded binary packages are in
> 
>
> C:\Users\stator-guest\AppData\Local\Temp\RtmpoVx40V\downloaded_packages
> 
>> library(TailRank)
> 
> Error: package or namespace load failed for TailRank in loadNamespace(i, 
> c(lib.loc, .libPaths()), versionCheck = vI[[i]]):
> 
> there is no package called Biobase
> 
> In addition: Warning message:
> 
> package TailRank was built under R version 3.4.2
> 
> 
> 
> 
>[[alternative HTML version deleted]]
> 
> __
> This email has been scanned for the NSW Ministry of Health by the Websense 
> Hosted Email Security System.
> Emails and attachments are monitored to ensure compliance with the NSW 
> Ministry of health's Electronic Messaging Policy.
> __
> ___
> Disclaimer: This message is intended for the addressee named and may contain 
> confidential information.
> If you are not the intended recipient, please delete it and notify the sender.
> Views expressed in this message are those of the individual sender, and are 
> not necessarily the views of the NSW Ministry of Health.
> ___
> This email has been scanned for the NSW Ministry of Health by the Websense 
> Hosted Email Security System.
> Emails and attachments are monitored to ensure compliance with the NSW 
> Ministry of Health's Electronic Messaging Policy.
> ___
> __
> R-help@r-project.org<mailto: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.

David Winsemius
Alameda, CA, USA

'Any technology distinguishable from magic is insufficiently advanced.'   
-Gehm's Corollary to Clarke's Third Law

__
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] "prob" package alternative

[David Winsemius]

> On Nov 1, 2017, at 12:51 PM, Tiby Kantrowitz  wrote:
> 
> The prob package has been archived because it depends upon some other
> packages which have issues.
> 
> However, such projects as Introduction to Probability and Statistics in R
> depend upon it for learning. There are a few other resources that also use
> it.
> 
> Does anyone know of any workarounds?
> 
> Someone at stack exchange mentioned using R 2.9.

I'm not sure I would trust that person. They seem a bit uninformed.

> However, that broke my
> RStudio (WSOD) and the dependent packages still wouldn't install, anyway.

The latest version of pkg-prob at the Archive directory of CRAN indicates that 
it was last updated within this year. The DESCRIPTION file indicates that it 
does not need compilation, but:

Depends: combinat, fAsianOptions

So there should be code in text files in its ../R directory which can be 
sourced from that directory.

~myuser_name$ ls /Users/../Downloads/prob/R 
characteristicfunctions.r   simulation.rutils-spaces.r
genData.R   spaces-examples.r   utils-subsets.r
misc.r  spaces-prob.r
prob.r  utils-events.r


Or you can install from source after downloading:

install.packages("~/Downloads/prob", repo=NULL,type="source")

# Success


>  library(prob)   # So does require having several other packages
Loading required package: combinat

Attaching package: ‘combinat’

The following object is masked from ‘package:utils’:

combn

Loading required package: fAsianOptions
Loading required package: timeDate

Attaching package: ‘timeDate’

The following object is masked from ‘package:cairoDevice’:

Cairo

The following objects are masked from ‘package:PerformanceAnalytics’:

kurtosis, skewness

Loading required package: timeSeries

Attaching package: ‘timeSeries’

The following object is masked from ‘package:zoo’:

time<-

Loading required package: fBasics


Rmetrics Package fBasics
Analysing Markets and calculating Basic Statistics
Copyright (C) 2005-2014 Rmetrics Association Zurich
Educational Software for Financial Engineering and Computational Science
Rmetrics is free software and comes with ABSOLUTELY NO WARRANTY.
https://www.rmetrics.org --- Mail to: i...@rmetrics.org
Loading required package: fOptions


Rmetrics Package fOptions
Pricing and Evaluating Basic Options
Copyright (C) 2005-2014 Rmetrics Association Zurich
Educational Software for Financial Engineering and Computational Science
Rmetrics is free software and comes with ABSOLUTELY NO WARRANTY.
https://www.rmetrics.org --- Mail to: i...@rmetrics.org

Attaching package: ‘prob’

The following objects are masked from ‘package:dplyr’:

intersect, setdiff, union

The following objects are masked from ‘package:base’:

intersect, setdiff, union

> 
> 
> 
> Tiby
> 
>   [[alternative HTML version deleted]]

A specific suggestion would be that you read the listinfo and the Posting Guide 
and learn to post in plain text.
> 
> __
> 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.

David Winsemius
Alameda, CA, USA

'Any technology distinguishable from magic is insufficiently advanced.'   
-Gehm's Corollary to Clarke's Third Law

__
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] ggplot inside function doesn't plot

[David Winsemius]

> On Nov 2, 2017, at 9:27 AM, Ed Siefker  wrote:
> 
> I have a function:
> 
> myplot <- function (X) {
>d <- plotCounts(dds2, gene=X, intgroup="condition", returnData=TRUE)
>png(paste("img/", X, ".png", sep=""))
>ggplot(d, aes(x=condition, y=count, color=condition)) +
>geom_point(position=position_jitter(w=0.1,h=0)) +
>scale_y_log10(breaks=c(25,100,400)) +
>ggtitle(X) +
>theme(plot.title = element_text(hjust = 0.5))
> 
>dev.off()
>}
> 
> 'd' is a dataframe
> 
> count  condition
> E11.5 F20HET BA40_quant   955.9788   E11.5 F20HET
> E11.5 F20HET BA45_quant   796.2863   E11.5 F20HET
> E11.5 F20HET BB84_quant   745.0340   E11.5 F20HET
> E11.5 F9.20DKO YEH3_quant 334.2994 E11.5 F9.20DKO
> E11.5 F9.20DKO fkm1_quant 313.7307 E11.5 F9.20DKO
> E11.5 F9.20DKO zzE2_quant 349.3313 E11.5 F9.20DKO
> 
> If I set X="Etv5" and paste the contents of the function into R, I get
> 'img/Etv5.png'
> If I run myplot(X), I get nothing.
> 
> 
>> X
> [1] "Etv5"
>> list.files("img")
> character(0)
>> myplot(X)
> null device
>  1
>> list.files("img")
> character(0)
>> d <- plotCounts(dds2, gene=X, intgroup="condition", returnData=TRUE)
>> png(paste("img/", X, ".png", sep=""))
>> ggplot(d, aes(x=condition, y=count, color=condition)) +
> + geom_point(position=position_jitter(w=0.1,h=0)) +
> + scale_y_log10(breaks=c(25,100,400)) +
> + ggtitle(X) +
> + theme(plot.title = element_text(hjust = 0.5))
>> dev.off()
> null device
>  1
>> list.files("img")
> [1] "Etv5.png"
> 
> Why doesn't my function work?

`ggplot` creates an object. You need to print it when used inside a function. 
Inside a function (in a more restricted environment) there is no 
parse-eval-print-loop.


> 
> __
> 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.

David Winsemius
Alameda, CA, USA

'Any technology distinguishable from magic is insufficiently advanced.'   
-Gehm's Corollary to Clarke's Third Law

__
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] ggplot inside function doesn't plot

[David Winsemius]

> On Nov 2, 2017, at 10:03 AM, Ed Siefker  wrote:
> 
> I don't really understand. I mean, I understand the solution is
> print(ggplot(...)).  But why is that required in a function and not at
> the console?

The REPL design of the interactive console is offered the user as a 
convenience, but I agree it's somewhat at odds with pure functional principles. 

> 
> Shouldn't I be able to rely on what I do at the console working in a
> script?  Is this inconsistent behavior by design?

By design? Yes. Completely by design. Functions are a method of encapsulating 
intermediate results in a manner that does not interact with the objects that 
exist outside the function. There are three plotting paradigms: base-graphics, 
lattice, and ggplot. The latter two are built with grid graphics and a both 
more functional and object oriented (loosely spoken since they return objects) 
than base graphics. Base graphics behaves the way you expect. Execution of 
`lines` or `points` will give you "pen-on-ink" output, actually pixel on device 
output.


> 
> 
> 
> On Thu, Nov 2, 2017 at 11:54 AM, David Winsemius  
> wrote:
>> 
>>> On Nov 2, 2017, at 9:27 AM, Ed Siefker  wrote:
>>> 
>>> I have a function:
>>> 
>>> myplot <- function (X) {
>>>   d <- plotCounts(dds2, gene=X, intgroup="condition", returnData=TRUE)
>>>   png(paste("img/", X, ".png", sep=""))
>>>   ggplot(d, aes(x=condition, y=count, color=condition)) +
>>>   geom_point(position=position_jitter(w=0.1,h=0)) +
>>>   scale_y_log10(breaks=c(25,100,400)) +
>>>   ggtitle(X) +
>>>   theme(plot.title = element_text(hjust = 0.5))
>>> 
>>>   dev.off()
>>>   }
>>> 
>>> 'd' is a dataframe
>>> 
>>>count  condition
>>> E11.5 F20HET BA40_quant   955.9788   E11.5 F20HET
>>> E11.5 F20HET BA45_quant   796.2863   E11.5 F20HET
>>> E11.5 F20HET BB84_quant   745.0340   E11.5 F20HET
>>> E11.5 F9.20DKO YEH3_quant 334.2994 E11.5 F9.20DKO
>>> E11.5 F9.20DKO fkm1_quant 313.7307 E11.5 F9.20DKO
>>> E11.5 F9.20DKO zzE2_quant 349.3313 E11.5 F9.20DKO
>>> 
>>> If I set X="Etv5" and paste the contents of the function into R, I get
>>> 'img/Etv5.png'
>>> If I run myplot(X), I get nothing.
>>> 
>>> 
>>>> X
>>> [1] "Etv5"
>>>> list.files("img")
>>> character(0)
>>>> myplot(X)
>>> null device
>>> 1
>>>> list.files("img")
>>> character(0)
>>>> d <- plotCounts(dds2, gene=X, intgroup="condition", returnData=TRUE)
>>>> png(paste("img/", X, ".png", sep=""))
>>>> ggplot(d, aes(x=condition, y=count, color=condition)) +
>>> + geom_point(position=position_jitter(w=0.1,h=0)) +
>>> + scale_y_log10(breaks=c(25,100,400)) +
>>> + ggtitle(X) +
>>> + theme(plot.title = element_text(hjust = 0.5))
>>>> dev.off()
>>> null device
>>> 1
>>>> list.files("img")
>>> [1] "Etv5.png"
>>> 
>>> Why doesn't my function work?
>> 
>> `ggplot` creates an object. You need to print it when used inside a 
>> function. Inside a function (in a more restricted environment) there is no 
>> parse-eval-print-loop.
>> 
>> 
>>> 
>>> ______
>>> 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.
>> 
>> David Winsemius
>> Alameda, CA, USA
>> 
>> 'Any technology distinguishable from magic is insufficiently advanced.'   
>> -Gehm's Corollary to Clarke's Third Law
>> 
>> 
>> 
>> 
>> 
> 
> __
> 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.

David Winsemius
Alameda, CA, USA

'Any technology distinguishable from magic is insufficiently advanced.'   
-Gehm's Corollary to Clarke's Third Law

__
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] "prob" package alternative

[David Winsemius]

> On Nov 2, 2017, at 11:15 AM, Tiby Kantrowitz  wrote:
> 
> The issue is fAsianOptions. Is there a version that works with the latest 
> version of R? If not, which version of it works with which version of R and 
> where can it be found? I tried several at the archive already.

sessionInfo()
R version 3.4.2 Patched (2017-10-04 r73465)
Platform: x86_64-apple-darwin15.6.0 (64-bit)
Running under: OS X El Capitan 10.11.6

> 


packageDescription("fAsianOptions")
Package: fAsianOptions
Version: 3010.79
Revision: 5522
Date: 2013-06-23
Title: EBM and Asian Option Valuation
Author: Diethelm Wuertz and many others, see the SOURCE file
Depends: R (>= 2.4.0), timeDate, timeSeries, fBasics, fOptions
Suggests: RUnit
Maintainer: Yohan Chalabi 
Description: Environment for teaching "Financial Engineering and Computational 
Finance"
Note: Several parts are still preliminary and may be changed in the future. this
  typically includes function and argument names, as well as defaults 
for
  arguments and return values.
LazyData: yes
License: GPL (>= 2)
URL: http://www.rmetrics.org
Packaged: 2013-06-23 18:22:14 UTC; yohan
NeedsCompilation: yes
Repository: CRAN
Date/Publication: 2013-06-24 01:53:27
Built: R 3.4.2; x86_64-apple-darwin15.6.0; 2017-11-01 22:45:02 UTC; unix



> Alternatively, is there another package that behaves similarly to prob?
> 
> 
> 
> On Wed, Nov 1, 2017 at 6:17 PM, David Winsemius  
> wrote:
> 
> > On Nov 1, 2017, at 12:51 PM, Tiby Kantrowitz  wrote:
> >
> > The prob package has been archived because it depends upon some other
> > packages which have issues.
> >
> > However, such projects as Introduction to Probability and Statistics in R
> > depend upon it for learning. There are a few other resources that also use
> > it.
> >
> > Does anyone know of any workarounds?
> >
> > Someone at stack exchange mentioned using R 2.9.
> 
> I'm not sure I would trust that person. They seem a bit uninformed.
> 
> > However, that broke my
> > RStudio (WSOD) and the dependent packages still wouldn't install, anyway.
> 
> The latest version of pkg-prob at the Archive directory of CRAN indicates 
> that it was last updated within this year. The DESCRIPTION file indicates 
> that it does not need compilation, but:
> 
> Depends: combinat, fAsianOptions
> 
> So there should be code in text files in its ../R directory which can be 
> sourced from that directory.
> 
> ~myuser_name$ ls /Users/../Downloads/prob/R
> characteristicfunctions.r   simulation.rutils-spaces.r
> genData.R   spaces-examples.r   
> utils-subsets.r
> misc.r  spaces-prob.r
> prob.r  utils-events.r
> 
> 
> Or you can install from source after downloading:
> 
> install.packages("~/Downloads/prob", repo=NULL,type="source")
> 
> # Success
> 
> 
> >  library(prob)   # So does require having several other packages
> Loading required package: combinat
> 
> Attaching package: ‘combinat’
> 
> The following object is masked from ‘package:utils’:
> 
> combn
> 
> Loading required package: fAsianOptions
> Loading required package: timeDate
> 
> Attaching package: ‘timeDate’
> 
> The following object is masked from ‘package:cairoDevice’:
> 
> Cairo
> 
> The following objects are masked from ‘package:PerformanceAnalytics’:
> 
> kurtosis, skewness
> 
> Loading required package: timeSeries
> 
> Attaching package: ‘timeSeries’
> 
> The following object is masked from ‘package:zoo’:
> 
> time<-
> 
> Loading required package: fBasics
> 
> 
> Rmetrics Package fBasics
> Analysing Markets and calculating Basic Statistics
> Copyright (C) 2005-2014 Rmetrics Association Zurich
> Educational Software for Financial Engineering and Computational Science
> Rmetrics is free software and comes with ABSOLUTELY NO WARRANTY.
> https://www.rmetrics.org --- Mail to: i...@rmetrics.org
> Loading required package: fOptions
> 
> 
> Rmetrics Package fOptions
> Pricing and Evaluating Basic Options
> Copyright (C) 2005-2014 Rmetrics Association Zurich
> Educational Software for Financial Engineering and Computational Science
> Rmetrics is free software and comes with ABSOLUTELY NO WARRANTY.
> https://www.rmetrics.org --- Mail to: i...@rmetrics.org
> 
> Attaching package: ‘prob’
> 
> The following objects are masked from ‘package:dplyr’:
> 
> intersect, setdiff, union
> 
> The following objects are masked from ‘package:base’:
> 
> intersect, setdiff, union
> 
> >
> >
> >
> > Tiby
> &g

Re: [R] "prob" package alternative

[David Winsemius]

> On Nov 2, 2017, at 12:07 PM, Tiby Kantrowitz  wrote:
> 
> Yes. That's the version I've been discussing that has non-zero exit status. 
> That situation is why CRAN retired the prob package. It's possible you 
> installed that library earlier in development and it's been "carried" along. 
> It no longer installs, now.
> 
> The problems with all of this seem to have started this month according to 
> the conversations. However, no one has mentioned any solutions or workarounds 
> except the one mentioned in passing (2.9).

Not true. Not even close to being true. I explained it all on the SO question 
page that you posted 2 days ago.

-- 

David.
> 
> Is there some other package that does something similar to prob that can be 
> used instead?
> 
> On Thu, Nov 2, 2017 at 2:29 PM, David Winsemius  
> wrote:
> 
> > On Nov 2, 2017, at 11:15 AM, Tiby Kantrowitz  wrote:
> >
> > The issue is fAsianOptions. Is there a version that works with the latest 
> > version of R? If not, which version of it works with which version of R and 
> > where can it be found? I tried several at the archive already.
> 
> sessionInfo()
> R version 3.4.2 Patched (2017-10-04 r73465)
> Platform: x86_64-apple-darwin15.6.0 (64-bit)
> Running under: OS X El Capitan 10.11.6
> 
> >
> 
> 
> packageDescription("fAsianOptions")
> Package: fAsianOptions
> Version: 3010.79
> Revision: 5522
> Date: 2013-06-23
> Title: EBM and Asian Option Valuation
> Author: Diethelm Wuertz and many others, see the SOURCE file
> Depends: R (>= 2.4.0), timeDate, timeSeries, fBasics, fOptions
> Suggests: RUnit
> Maintainer: Yohan Chalabi 
> Description: Environment for teaching "Financial Engineering and 
> Computational Finance"
> Note: Several parts are still preliminary and may be changed in the future. 
> this
>   typically includes function and argument names, as well as defaults 
> for
>   arguments and return values.
> LazyData: yes
> License: GPL (>= 2)
> URL: http://www.rmetrics.org
> Packaged: 2013-06-23 18:22:14 UTC; yohan
> NeedsCompilation: yes
> Repository: CRAN
> Date/Publication: 2013-06-24 01:53:27
> Built: R 3.4.2; x86_64-apple-darwin15.6.0; 2017-11-01 22:45:02 UTC; unix
> 
> 
> 
> > Alternatively, is there another package that behaves similarly to prob?
> >
> >
> >
> > On Wed, Nov 1, 2017 at 6:17 PM, David Winsemius  
> > wrote:
> >
> > > On Nov 1, 2017, at 12:51 PM, Tiby Kantrowitz  wrote:
> > >
> > > The prob package has been archived because it depends upon some other
> > > packages which have issues.
> > >
> > > However, such projects as Introduction to Probability and Statistics in R
> > > depend upon it for learning. There are a few other resources that also use
> > > it.
> > >
> > > Does anyone know of any workarounds?
> > >
> > > Someone at stack exchange mentioned using R 2.9.
> >
> > I'm not sure I would trust that person. They seem a bit uninformed.
> >
> > > However, that broke my
> > > RStudio (WSOD) and the dependent packages still wouldn't install, anyway.
> >
> > The latest version of pkg-prob at the Archive directory of CRAN indicates 
> > that it was last updated within this year. The DESCRIPTION file indicates 
> > that it does not need compilation, but:
> >
> > Depends: combinat, fAsianOptions
> >
> > So there should be code in text files in its ../R directory which can be 
> > sourced from that directory.
> >
> > ~myuser_name$ ls /Users/../Downloads/prob/R
> > characteristicfunctions.r   simulation.r
> > utils-spaces.r
> > genData.R   spaces-examples.r   
> > utils-subsets.r
> > misc.r  spaces-prob.r
> > prob.r  utils-events.r
> >
> >
> > Or you can install from source after downloading:
> >
> > install.packages("~/Downloads/prob", repo=NULL,type="source")
> >
> > # Success
> >
> >
> > >  library(prob)   # So does require having several other packages
> > Loading required package: combinat
> >
> > Attaching package: ‘combinat’
> >
> > The following object is masked from ‘package:utils’:
> >
> > combn
> >
> > Loading required package: fAsianOptions
> > Loading required package: timeDate
> >
> > Attaching package: ‘timeDate’
> >
> > The following object is masked from ‘package:cairoDevice’:
> >

Re: [R] "prob" package alternative

[David Winsemius]

> On Nov 2, 2017, at 1:09 PM, Tiby Kantrowitz  wrote:
> 
> Yes, that is exactly what I was doing two days ago.
> 
> Warning in install.packages :
>   installation of package ‘fAsianOptions_3010.79.tar.gz’ had non-zero exit 
> status
> 
> Which is what a reading of the explanation for why "prob" was retired leads 
> one to expect. Do you have some other suggestion about how to get it to work? 
> I notice you're not using Windows which might have a relationship to why it 
> is working for you.


I explained what you needed to do if you were on Windows. You do need to read 
the explanation more closely, identify the parts you don't understand (in all 
likelihood the point I made about needing the proper development tools for 
Windows) _and_ read:

https://cran.r-project.org/doc/manuals/r-release/R-admin.html#Installing-packages

If you have trouble understanding or getting success after thoroughly reading 
these directions, you need to explain what you have done (preserving the order 
of all operations) and post a complete transcript of all commands and complete 
error messages.



-- 
David.
> Otherwise, do you know of some other package that could be used instead of 
> prob?
> 
> On Thu, Nov 2, 2017 at 3:44 PM, David Winsemius  
> wrote:
> 
> > On Nov 2, 2017, at 12:07 PM, Tiby Kantrowitz  wrote:
> >
> > Yes. That's the version I've been discussing that has non-zero exit status. 
> > That situation is why CRAN retired the prob package. It's possible you 
> > installed that library earlier in development and it's been "carried" 
> > along. It no longer installs, now.
> >
> > The problems with all of this seem to have started this month according to 
> > the conversations. However, no one has mentioned any solutions or 
> > workarounds except the one mentioned in passing (2.9).
> 
> Not true. Not even close to being true. I explained it all on the SO question 
> page that you posted 2 days ago.
> 
> --
> 
> David.
> >
> > Is there some other package that does something similar to prob that can be 
> > used instead?
> >
> > On Thu, Nov 2, 2017 at 2:29 PM, David Winsemius  
> > wrote:
> >
> > > On Nov 2, 2017, at 11:15 AM, Tiby Kantrowitz  wrote:
> > >
> > > The issue is fAsianOptions. Is there a version that works with the latest 
> > > version of R? If not, which version of it works with which version of R 
> > > and where can it be found? I tried several at the archive already.
> >
> > sessionInfo()
> > R version 3.4.2 Patched (2017-10-04 r73465)
> > Platform: x86_64-apple-darwin15.6.0 (64-bit)
> > Running under: OS X El Capitan 10.11.6
> >
> > >
> >
> >
> > packageDescription("fAsianOptions")
> > Package: fAsianOptions
> > Version: 3010.79
> > Revision: 5522
> > Date: 2013-06-23
> > Title: EBM and Asian Option Valuation
> > Author: Diethelm Wuertz and many others, see the SOURCE file
> > Depends: R (>= 2.4.0), timeDate, timeSeries, fBasics, fOptions
> > Suggests: RUnit
> > Maintainer: Yohan Chalabi 
> > Description: Environment for teaching "Financial Engineering and 
> > Computational Finance"
> > Note: Several parts are still preliminary and may be changed in the future. 
> > this
> >   typically includes function and argument names, as well as 
> > defaults for
> >   arguments and return values.
> > LazyData: yes
> > License: GPL (>= 2)
> > URL: http://www.rmetrics.org
> > Packaged: 2013-06-23 18:22:14 UTC; yohan
> > NeedsCompilation: yes
> > Repository: CRAN
> > Date/Publication: 2013-06-24 01:53:27
> > Built: R 3.4.2; x86_64-apple-darwin15.6.0; 2017-11-01 22:45:02 UTC; unix
> >
> >
> >
> > > Alternatively, is there another package that behaves similarly to prob?
> > >
> > >
> > >
> > > On Wed, Nov 1, 2017 at 6:17 PM, David Winsemius  
> > > wrote:
> > >
> > > > On Nov 1, 2017, at 12:51 PM, Tiby Kantrowitz  wrote:
> > > >
> > > > The prob package has been archived because it depends upon some other
> > > > packages which have issues.
> > > >
> > > > However, such projects as Introduction to Probability and Statistics in 
> > > > R
> > > > depend upon it for learning. There are a few other resources that also 
> > > > use
> > > > it.
> > > >
> > > > Does anyone know of any workarounds?
> > > >
> > > > Someone at stack exchange mentioned using R 2.

Re: [R] "prob" package alternative

[David Winsemius]

> On Nov 2, 2017, at 2:14 PM, Tiby Kantrowitz  wrote:
> 
> Rtools is not available for the current version of R.

Really? If true, I'm surprised and not able to help. I do see an Rtools34.exe 
at https://cran.r-project.org/bin/windows/Rtools/
-- 
David.

> 
> What I'm looking for is an alternative package or how others have managed to 
> create workarounds.
> 
> On Thu, Nov 2, 2017 at 4:25 PM, David Winsemius  
> wrote:
> 
> > On Nov 2, 2017, at 1:09 PM, Tiby Kantrowitz  wrote:
> >
> > Yes, that is exactly what I was doing two days ago.
> >
> > Warning in install.packages :
> >   installation of package ‘fAsianOptions_3010.79.tar.gz’ had non-zero exit 
> > status
> >
> > Which is what a reading of the explanation for why "prob" was retired leads 
> > one to expect. Do you have some other suggestion about how to get it to 
> > work? I notice you're not using Windows which might have a relationship to 
> > why it is working for you.
> 
> 
> I explained what you needed to do if you were on Windows. You do need to read 
> the explanation more closely, identify the parts you don't understand (in all 
> likelihood the point I made about needing the proper development tools for 
> Windows) _and_ read:
> 
> https://cran.r-project.org/doc/manuals/r-release/R-admin.html#Installing-packages
> 
> If you have trouble understanding or getting success after thoroughly reading 
> these directions, you need to explain what you have done (preserving the 
> order of all operations) and post a complete transcript of all commands and 
> complete error messages.
> 
> 
> 
> --
> David.
> > Otherwise, do you know of some other package that could be used instead of 
> > prob?
> >
> > On Thu, Nov 2, 2017 at 3:44 PM, David Winsemius  
> > wrote:
> >
> > > On Nov 2, 2017, at 12:07 PM, Tiby Kantrowitz  wrote:
> > >
> > > Yes. That's the version I've been discussing that has non-zero exit 
> > > status. That situation is why CRAN retired the prob package. It's 
> > > possible you installed that library earlier in development and it's been 
> > > "carried" along. It no longer installs, now.
> > >
> > > The problems with all of this seem to have started this month according 
> > > to the conversations. However, no one has mentioned any solutions or 
> > > workarounds except the one mentioned in passing (2.9).
> >
> > Not true. Not even close to being true. I explained it all on the SO 
> > question page that you posted 2 days ago.
> >
> > --
> >
> > David.
> > >
> > > Is there some other package that does something similar to prob that can 
> > > be used instead?
> > >
> > > On Thu, Nov 2, 2017 at 2:29 PM, David Winsemius  
> > > wrote:
> > >
> > > > On Nov 2, 2017, at 11:15 AM, Tiby Kantrowitz  wrote:
> > > >
> > > > The issue is fAsianOptions. Is there a version that works with the 
> > > > latest version of R? If not, which version of it works with which 
> > > > version of R and where can it be found? I tried several at the archive 
> > > > already.
> > >
> > > sessionInfo()
> > > R version 3.4.2 Patched (2017-10-04 r73465)
> > > Platform: x86_64-apple-darwin15.6.0 (64-bit)
> > > Running under: OS X El Capitan 10.11.6
> > >
> > > >
> > >
> > >
> > > packageDescription("fAsianOptions")
> > > Package: fAsianOptions
> > > Version: 3010.79
> > > Revision: 5522
> > > Date: 2013-06-23
> > > Title: EBM and Asian Option Valuation
> > > Author: Diethelm Wuertz and many others, see the SOURCE file
> > > Depends: R (>= 2.4.0), timeDate, timeSeries, fBasics, fOptions
> > > Suggests: RUnit
> > > Maintainer: Yohan Chalabi 
> > > Description: Environment for teaching "Financial Engineering and 
> > > Computational Finance"
> > > Note: Several parts are still preliminary and may be changed in the 
> > > future. this
> > >   typically includes function and argument names, as well as 
> > > defaults for
> > >   arguments and return values.
> > > LazyData: yes
> > > License: GPL (>= 2)
> > > URL: http://www.rmetrics.org
> > > Packaged: 2013-06-23 18:22:14 UTC; yohan
> > > NeedsCompilation: yes
> > > Repository: CRAN
> > > Date/Publication: 2013-06-24 01:53:27
> > > Built: R 3.4.2; x86_64-apple-dar

Re: [R] "prob" package alternative

[David Winsemius]

> On Nov 2, 2017, at 3:46 PM, Tiby Kantrowitz  wrote:
> 
> Thanks. I found that, and installed it and got the same message. Here:
> 
> RTools version 3.4 
> 
> install.packages("fAsianOptions_3010.tar.gz",

I don't see a path to that file's location.

The expansion from pkg_version.tar.gz might be automatic, but I do generally 
expand such file into its own directory and then drag that directory icon into 
the console so my command would look like:

install.packages("~/Downloads/fAsianOptions", dependencies=TRUE, repos=NULL, 
type = "source") 


I suggest you learn to post in plain text since your HTML postings are loosing 
their carriage returns.

-- 
David




> dependencies=TRUE, repos=NULL, type = "source") 
> Installing package into ‘C:/Users/Tlk7/Documents/R/win-library/3.4’ (as ‘lib’ 
> is unspecified) Warning: invalid package 'fAsianOptions_3010.tar.gz' Error: 
> ERROR: no packages specified Warning in install.packages : running command 
> '"C:/PROGRA~1/R/R-34~1.2/bin/x64/R" CMD INSTALL -l 
> "C:\Users\Tlk7\Documents\R\win-library\3.4" "fAsianOptions_3010.tar.gz"' had 
> status 1 Warning in install.packages : installation of package 
> ‘fAsianOptions_3010.tar.gz’ had non-zero exit status
> 
> On Thu, Nov 2, 2017 at 6:06 PM, David Winsemius  
> wrote:
> 
> > On Nov 2, 2017, at 2:14 PM, Tiby Kantrowitz  wrote:
> >
> > Rtools is not available for the current version of R.
> 
> Really? If true, I'm surprised and not able to help. I do see an Rtools34.exe 
> at https://cran.r-project.org/bin/windows/Rtools/
> --
> David.
> 
> >
> > What I'm looking for is an alternative package or how others have managed 
> > to create workarounds.
> >
> > On Thu, Nov 2, 2017 at 4:25 PM, David Winsemius  
> > wrote:
> >
> > > On Nov 2, 2017, at 1:09 PM, Tiby Kantrowitz  wrote:
> > >
> > > Yes, that is exactly what I was doing two days ago.
> > >
> > > Warning in install.packages :
> > >   installation of package ‘fAsianOptions_3010.79.tar.gz’ had non-zero 
> > > exit status
> > >
> > > Which is what a reading of the explanation for why "prob" was retired 
> > > leads one to expect. Do you have some other suggestion about how to get 
> > > it to work? I notice you're not using Windows which might have a 
> > > relationship to why it is working for you.
> >
> >
> > I explained what you needed to do if you were on Windows. You do need to 
> > read the explanation more closely, identify the parts you don't understand 
> > (in all likelihood the point I made about needing the proper development 
> > tools for Windows) _and_ read:
> >
> > https://cran.r-project.org/doc/manuals/r-release/R-admin.html#Installing-packages
> >
> > If you have trouble understanding or getting success after thoroughly 
> > reading these directions, you need to explain what you have done 
> > (preserving the order of all operations) and post a complete transcript of 
> > all commands and complete error messages.
> >
> >
> >
> > --
> > David.
> > > Otherwise, do you know of some other package that could be used instead 
> > > of prob?
> > >
> > > On Thu, Nov 2, 2017 at 3:44 PM, David Winsemius  
> > > wrote:
> > >
> > > > On Nov 2, 2017, at 12:07 PM, Tiby Kantrowitz  wrote:
> > > >
> > > > Yes. That's the version I've been discussing that has non-zero exit 
> > > > status. That situation is why CRAN retired the prob package. It's 
> > > > possible you installed that library earlier in development and it's 
> > > > been "carried" along. It no longer installs, now.
> > > >
> > > > The problems with all of this seem to have started this month according 
> > > > to the conversations. However, no one has mentioned any solutions or 
> > > > workarounds except the one mentioned in passing (2.9).
> > >
> > > Not true. Not even close to being true. I explained it all on the SO 
> > > question page that you posted 2 days ago.
> > >
> > > --
> > >
> > > David.
> > > >
> > > > Is there some other package that does something similar to prob that 
> > > > can be used instead?
> > > >
> > > > On Thu, Nov 2, 2017 at 2:29 PM, David Winsemius 
> > > >  wrote:
> > > >
> > > > > On Nov 2, 2017, at 11:15 AM, Tiby Kantrowitz  
> >

Re: [R] Problem with r project in ubuntu xenial

[David Winsemius]

> On Nov 3, 2017, at 5:09 PM, peter dalgaard  wrote:
> 
> 
>> On 3 Nov 2017, at 23:39 , George Balas  wrote:
>> 
>> I have a problem with R in Ubuntu 16.04. I do not know if it is mine pc or
>> general problem but I was not able to find solution on Internet.
>> First of all I can not change locale to greek by getting this message:
>> "In Sys.setlocale("LC_CTYPE", "Greek") :
>> OS reports request to set locale to "Greek" cannot be honored"
> 
> The Greek locale is likely not called "Greek" outside of Windows. More likely 
> "el_GR.UTF-8" or thereabouts (check your locale database, I'm on a Mac). 
> These things are not standardized across platforms.
> 

Also

>From help(locales): Attempts to change the character set by 
>Sys.setlocale("LC_CTYPE") that implies a different character set during a 
>session may not work and are likely to lead to some confusion because it may 
>not affect the native encoding.


>> Second and more serious is that I can not use some functions like
>> graph_from_adjacency_matrix or print_all I get these messeges:
>> "could not find function "graph_from_adjacency_matrix""
>> "could not find function "print_all"".
> 
> Missing library(igraph)?
> 
> -pd
> 
>> I am using R version 3.4.2 (2017-09-28) -- "Short Summer" either on rstudio
>> or ubuntu terminal.
>> On my pc I also run win 10 with the same installs and I do not have the
>> above problems, but I work on ubuntu and can not change Os all the time.
>> Please help me.
>> 
>> Thank you for your time,
>> George
>> gbala...@gmail.com
>> 
>>  [[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.
> 
> -- 
> Peter Dalgaard, Professor,
> Center for Statistics, Copenhagen Business School
> Solbjerg Plads 3, 2000 Frederiksberg, Denmark
> Phone: (+45)38153501
> Office: A 4.23
> Email: pd@cbs.dk  Priv: pda...@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.

David Winsemius
Alameda, CA, USA

'Any technology distinguishable from magic is insufficiently advanced.'   
-Gehm's Corollary to Clarke's Third Law

__
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] Change colour of line in logi.hist.plot

[David Winsemius]

> On Nov 5, 2017, at 6:03 AM, Göran Bergqvist 
>  wrote:
> 
> I am using the function logi.hist.plot in package popbio. I want to change 
> the colour of the probability line from the default red to black. I have not 
> been able to find out how to do that.

If you look at the code for that function you see this line:

logi.curve <- function(independ, depend, mod = logi.mod, 
col.cur = "red", lwd.cur = 4) {
Since it's the only occurence of "red" that's probably where "the money lies". 
You can either do a "hard hack" where you alter the value of the parameters to 
that inner function,

... or you can add a named parameter after the dots in the outer parameter list 
such as:
   

 function (independ, depend, logi.mod = 1, type = "dit", boxp = TRUE, 
rug = FALSE, ylabel = "Probability", ylabel2 = "Frequency", 
xlabel = "", mainlabel = "", las.h = 1, counts = FALSE, ...,  col.cur = 
"red")

... and modify the inner function to read:

 logi.curve <- function(independ, depend, mod = logi.mod, 
col.cur = col.cur, lwd.cur = 4) {


>   [[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.

David Winsemius
Alameda, CA, USA

'Any technology distinguishable from magic is insufficiently advanced.'   
-Gehm's Corollary to Clarke's Third Law

__
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] Help in R

[David Winsemius]

> On Nov 5, 2017, at 9:28 AM, Ahsan Zahir via R-help  
> wrote:
> 
> 
> Hey,
> 
> I am a beginner in R.
> 
> How do I read last 10 values from column- Movie, from a dataset?

Some questions are so simple that they strongly suggest no prior effort at 
self-leanrning. In such cases the usual recommendation given at Rhelp is that 
you read an introductory text. Many of us used the "Introduction to R" that is 
shipped with every copy of R:

https://cran.r-project.org/doc/manuals/r-release/R-intro.pdf


> 
> Pls help.
> 
> Sent from my iPhone
> 
> __
> 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.

David Winsemius
Alameda, CA, USA

'Any technology distinguishable from magic is insufficiently advanced.'   
-Gehm's Corollary to Clarke's Third Law

__
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] Survival model error

[David Winsemius]

You should stop trying to use matrices on the RHS and using separate vectors to 
Surv. Instead use a data argument and have the names in your formula refer to 
column names. 

— 
David

Sent from my iPhone

> On Nov 5, 2017, at 7:21 PM, Meghna Govil via R-help  
> wrote:
> 
> 
> Hi - Below is my code and then the error when I run the last line. 
> 
> time_np <- train1_na$tte
> event_np <- train1_na$censored
> 
> 
> X_np <- cbind(
>  train1_na$AMT, 
>  train1_na$DISCOUNT_AMT,
>  train1_na$high_price_pcnt,
>  train1_na$EM_RECEIVED,
>  train1_na$DM_RECEIVED,
>  train1_na$TXN_WITH_RINGCODE,
>  train1_na$WEB,
>  train1_na$clearance_pcnt,
>  train1_na$bts_pcnt,
>  train1_na$sales_pcnt,
>  train1_na$holiday_pcnt,
>  train1_na$TXN,
>  train1_na$REDEEMED_REWARDS
> 
>  )
> 
> # Kaplan-Meier non-parametric analysis 
> 
> kmsurvival_np <- survfit(Surv(time_np,event_np) ~ X_np)
> 
> Error in `[.default`(y, who, 1) : (subscript) logical subscript too long
> 
> Any ideas?  
> I have tried several things and still get this error.  
> 
> 
> Thanks,
> Meghna
> 
> __
> 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] Survival model error

[David Winsemius]

> On Nov 6, 2017, at 5:45 AM, Meghna Govil  wrote:
> 
> Thanks David. Could you show me how to do that in my example ? 

Possibly:

kmsurvival_np <- survfit(Surv( tte, censored) ~ . , data=train1_na)

I say "possibly" because I don't know whether all the columns of `train1_na` 
were included in the matrix you constructed. If they weren't then you would 
need to do something like:

kmsurvival_np <- survfit(Surv( tte, censored) ~  , data=train1_na)


The formalism of `reg_func( Y , Xm)` Where Xm is a matrix is typical for 
various machine-learning types of procedures, but breaks the 
"process-names-in-an-environment model of typical R regression functions. R 
regression functions often omit "step-down" options so loved by beginning SAS 
users.


> 
> Thanks,
> Meghna
> 
>> On Nov 6, 2017, at 12:58 AM, David Winsemius  wrote:
>> 
>> You should stop trying to use matrices on the RHS and using separate vectors 
>> to Surv. Instead use a data argument and have the names in your formula 
>> refer to column names. 
>> 
>> — 
>> David
>> 
>> Sent from my iPhone
>> 
>>> On Nov 5, 2017, at 7:21 PM, Meghna Govil via R-help  
>>> wrote:
>>> 
>>> 
>>> Hi - Below is my code and then the error when I run the last line. 
>>> 
>>> time_np <- train1_na$tte
>>> event_np <- train1_na$censored
>>> 
>>> 
>>> X_np <- cbind(
>>> train1_na$AMT, 
>>> train1_na$DISCOUNT_AMT,
>>> train1_na$high_price_pcnt,
>>> train1_na$EM_RECEIVED,
>>> train1_na$DM_RECEIVED,
>>> train1_na$TXN_WITH_RINGCODE,
>>> train1_na$WEB,
>>> train1_na$clearance_pcnt,
>>> train1_na$bts_pcnt,
>>> train1_na$sales_pcnt,
>>> train1_na$holiday_pcnt,
>>> train1_na$TXN,
>>> train1_na$REDEEMED_REWARDS
>>> 
>>> )
>>> 
>>> # Kaplan-Meier non-parametric analysis 
>>> 
>>> kmsurvival_np <- survfit(Surv(time_np,event_np) ~ X_np)
>>> 
>>> Error in `[.default`(y, who, 1) : (subscript) logical subscript too long
>>> 
>>> Any ideas?  
>>> I have tried several things and still get this error.  
>>> 
>>> 
>>> Thanks,
>>> Meghna
>>> 
>>> __
>>> 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.
>> 
> 

David Winsemius
Alameda, CA, USA

'Any technology distinguishable from magic is insufficiently advanced.'   
-Gehm's Corollary to Clarke's Third Law

__
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] Fitdistrplus and Custom Probability Density

[David Winsemius]

> On Nov 7, 2017, at 6:58 AM, Eric Berger  wrote:
> 
> Why not define your own functions based on d?
> e.g.
> myCumDist <- function(x) { integrate(d, lower=-Inf, upper=x)$value  }
> myQuantile <- function(x) { uniroot(f=function(y) { h(y) - x },
> interval=c(-5,5)) }  # limits -5,5 should be replaced by your own which
> might require some fiddling

My test gave an error regarding a missing `h` function. Perhaps you meant 
something like:

myQuantile <- function(x) { uniroot(f=function(y) { myCumDist(y)-x },
   interval=c(-5,5))$root }

I added the $root extraction so it now returns a numeric value instead of a 
list. I wondered if that might be somewhat slow, and considered the possibility 
of building a close aproximation that didn't require a combined integration and 
rootfinding operation for every value using approxfun() on a suitable range of 
p-values in the range [0,1]

-- 
David.
> 
> e.g.
> d <- function(x) { exp(-x^2/2)/(sqrt(2*pi)) }  # just an example for you to
> test with; use your own density d(x) in your case
> 
> Then define myCumDist, myQuantile as above and compare with pnorm, qnorm.
> 
> HTH,
> Eric
> 
> 
> 
> 
> On Tue, Nov 7, 2017 at 4:22 PM, Lorenzo Isella 
> wrote:
> 
>> Dear All,
>> Apologies for not providing a reproducible example, but if I could, then I
>> would be able to answer myself my question.
>> Essentially, I am trying to fit a very complicated custom probability
>> distribution to some data.
>> Fitdistrplus does in principle everything which I need, but if require me
>> to specify not only the density function d, but also the cumulative p and
>> and inverse cumulative function q (see for instance
>> 
>> http://www.stat.umn.edu/geyer/old/5101/rlook.html
>> 
>> to understand what these quantities are in the case of a normal
>> distribution).
>> 
>> The analytical calculation of p and q is a big task in my case, so my
>> question is if there is a workaround for this, i.e. a way to fit the
>> unknown parameters of my probability distribution without specifying (at
>> least analytically) p and q, but only the density d.
>> Many thanks
>> 
>> Lorenzo
>> 
>>[[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.

David Winsemius
Alameda, CA, USA

'Any technology distinguishable from magic is insufficiently advanced.'   
-Gehm's Corollary to Clarke's Third Law

__
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] Problem with r project in ubuntu xenial

[David Winsemius]

> On Nov 7, 2017, at 3:46 PM, George Balas  wrote:
> 
> For anyone who sees this conversation.
> 
> There is a bug in installation of igraph in R language in Ubuntu. There is
> a solution in stackoverflow. We have to use the devtools. Write this code:
> install.packages("devtools")
> library(devtools)
> install_github("igraph/rigraph")

If there is a bug in the development version of igraph (which is not really 
part of R but rather a contributed package) then the correct place to send a 
message would have been, not to rhelp, but rather to the maintainer. It would 
of course be necessary to say what the bug appears to be. I don't see any error 
message or description of a bug in this message. The package DESCRIPTION files 
says:

BugReports: https://github.com/igraph/igraph/issues

The github page says the proper installation sequence for the development 
version is:

devtools::install_github("gaborcsardi/pkgconfig")
devtools::install_github("igraph/rigraph")

Rhelp is not the correct place for pre-release bug discussion of contributed 
packages.

Again, the github pages suggests that you send issues to the correct mailing 
list: igraph-help mailing list

https://lists.nongnu.org/mailman/listinfo/igraph-help


-- 
David.
> 
> If there are errors installing devtools just install any package that
> comments.
> 
> On Nov 5, 2017 00:07, "George Balas"  wrote:
> 
>> -Well it seems that it is getting "el_GR.UTF-8" but still I am not able
>> to read files written in greek, there are only "" instead of letters.
>> -Also, I forgot to mention that I do load igraph library when I try 
>> "graph_from_adjacency_matrix".
>> When I check igraph in packages dialog I can not see functions with
>> underscores between words, only dots.
>> 
>> 2017-11-04 2:22 GMT+02:00 David Winsemius :
>> 
>>> 
>>>> On Nov 3, 2017, at 5:09 PM, peter dalgaard  wrote:
>>>> 
>>>> 
>>>>> On 3 Nov 2017, at 23:39 , George Balas  wrote:
>>>>> 
>>>>> I have a problem with R in Ubuntu 16.04. I do not know if it is mine
>>> pc or
>>>>> general problem but I was not able to find solution on Internet.
>>>>> First of all I can not change locale to greek by getting this message:
>>>>> "In Sys.setlocale("LC_CTYPE", "Greek") :
>>>>> OS reports request to set locale to "Greek" cannot be honored"
>>>> 
>>>> The Greek locale is likely not called "Greek" outside of Windows. More
>>> likely "el_GR.UTF-8" or thereabouts (check your locale database, I'm on a
>>> Mac). These things are not standardized across platforms.
>>>> 
>>> 
>>> Also
>>> 
>>> From help(locales): Attempts to change the character set by
>>> Sys.setlocale("LC_CTYPE") that implies a different character set during a
>>> session may not work and are likely to lead to some confusion because it
>>> may not affect the native encoding.
>>> 
>>> 
>>>>> Second and more serious is that I can not use some functions like
>>>>> graph_from_adjacency_matrix or print_all I get these messeges:
>>>>> "could not find function "graph_from_adjacency_matrix""
>>>>> "could not find function "print_all"".
>>>> 
>>>> Missing library(igraph)?
>>>> 
>>>> -pd
>>>> 
>>>>> I am using R version 3.4.2 (2017-09-28) -- "Short Summer" either on
>>> rstudio
>>>>> or ubuntu terminal.
>>>>> On my pc I also run win 10 with the same installs and I do not have the
>>>>> above problems, but I work on ubuntu and can not change Os all the
>>> time.
>>>>> Please help me.
>>>>> 
>>>>> Thank you for your time,
>>>>> George
>>>>> gbala...@gmail.com
>>>>> 
>>>>> [[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/posti
>>> ng-guide.html
>>>>> and provide commented, minimal, self-contained, reproducible code.
>>>> 
>>>> --
>>>> Peter Dalgaard, Professor,
>>>> Center for St

Re: [R] Primer for working with survey data in R

[David Winsemius]

> On Nov 11, 2017, at 11:56 AM, Kevin Taylor  wrote:
> 
> I am taking a behavioral stats graduate class and the instructor is using
> SPSS. I'm trying to follow along in R.
> 
> Recently in class we started working with scales and survey data, computing
> Cronbach's Alpha, reversing values for reverse coded items, etc.
> 
> Also, SPSS has some built in functionality for entering the meta-data for
> your survey, e.g. the possible values for items, the text of the question,
> etc.
> 
> I haven't been able to find any survey guidance for R other than how to run
> the actual calculations (Cronbach's, reversing values).
> 
> Are there tutorials, books, or other primers, that would guide a newbie
> step by step through using R for working with survey data? It would be
> helpful to see how others are doing these things. (Not just how to run the
> mathematical operations but how to work with and manage the data.) Possibly
> this would be in conjunction with some packages such as Likert or Scales.

Try looking at:

http://personality-project.org/r/psych/

-- 
David.
> 
> TIA.
> 
> --Kevin
> 
>   [[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.

David Winsemius
Alameda, CA, USA

'Any technology distinguishable from magic is insufficiently advanced.'   
-Gehm's Corollary to Clarke's Third Law

__
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] Primer for working with survey data in R

[David Winsemius]

> On Nov 11, 2017, at 2:36 PM, Fox, John  wrote:
> 
> Dear Kevin,
> 
> In addition to the advice you've received, take a look at the survey package. 
> It's not quite what you're asking for, but in fact it's probably more useful, 
> in that it provides correct statistical inference for data collected in 
> complex surveys. The package is described in an article,  T. Lumley (2004), 
> Analysis of complex survey samples, Journal of Statistical Software 9(1): 
> 1-19, and a book, T. Lumley, Complex Surveys: A Guide to Analysis Using R, 
> Wiley, 2010, both by the package author.

Although the same thought occurred to me after reading the initial question, I 
decided against suggesting the survey package. I consulted the recommended book 
above and it has none of the requested statistics. It is designed for surveys 
that use complex sampling designs requiring weighting the observations. I also 
consulted the Social Sciences Task View and it seemed unhelpful for the 
specific requests. It seemed likely to me that even a graduate course in 
behavioral statistics would be focussed on the sorts of questions that the 
psych package delivers. The  website maintained by Revelle has several 
tutorials that include developed examples using R to deliver the requested 
measure. Obviously "reversing values" is something that would require learning 
basic R manipulation of factor variables.

-- 
David.

> 
> I hope that this helps,
> John
> 
>> -Original Message-
>> From: R-help [mailto:r-help-boun...@r-project.org] On Behalf Of Kevin Taylor
>> Sent: Saturday, November 11, 2017 2:57 PM
>> To: r-help@r-project.org
>> Subject: [R] Primer for working with survey data in R
>> 
>> I am taking a behavioral stats graduate class and the instructor is using 
>> SPSS.
>> I'm trying to follow along in R.
>> 
>> Recently in class we started working with scales and survey data, computing
>> Cronbach's Alpha, reversing values for reverse coded items, etc.
>> 
>> Also, SPSS has some built in functionality for entering the meta-data for 
>> your
>> survey, e.g. the possible values for items, the text of the question, etc.
>> 
>> I haven't been able to find any survey guidance for R other than how to run 
>> the
>> actual calculations (Cronbach's, reversing values).
>> 
>> Are there tutorials, books, or other primers, that would guide a newbie step 
>> by
>> step through using R for working with survey data? It would be helpful to see
>> how others are doing these things. (Not just how to run the mathematical
>> operations but how to work with and manage the data.) Possibly this would be
>> in conjunction with some packages such as Likert or Scales.
>> 
>> TIA.
>> 
>> --Kevin
>> 
>>  [[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.

David Winsemius
Alameda, CA, USA

'Any technology distinguishable from magic is insufficiently advanced.'   
-Gehm's Corollary to Clarke's Third Law

__
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] Bootstrap analysis from a conditional logistic regression

[David Winsemius]

> On Nov 13, 2017, at 2:01 PM, Nelly Reduan  wrote:
> 
> Nelly Reduan a partag� un fichier OneDrive avec vous. Pour l�afficher, 
> cliquez sur le lien ci-dessous.
> 
> 
> <https://1drv.ms/u/s!Apkg2VlgfYyDgRAeVIM0nEajx0Fb>
> [https://r1.res.office365.com/owa/prem/images/dc-png_20.png]<https://1drv.ms/u/s!Apkg2VlgfYyDgRAeVIM0nEajx0Fb>
> 
> Screenshot 2017-11-12 
> 18.49.43.png<https://1drv.ms/u/s!Apkg2VlgfYyDgRAeVIM0nEajx0Fb>
> 
> 
> 
> 
> Hello
> 
> How can I perform a bootstrap analysis from a conditional logistic 
> regression? The model has been built using the `clogit` function (`survival` 
> package)? The model has the following structure:
> 
>mod <- clogit(event ~ forest + log_area +forest:log_time  + 
> cluster(ID_individual)  +   strata(ID_strata), method = "efron", data = data 
> , x=T, y=T)
> 
> Using bootstrapping, I would like to have a measure of uncertainty around the 
> estimates of beta coefficients.
> 
> I am using the following code but I don't know how to consider strata and 
> cluster arguments.
> 
>library(boot)
>boot.clogit <- function(data, indices){
>  new_data <- data[indices,]
>  mod <- clogit(event ~ forest + log_area + forest:log_time  + 
> cluster(ID_individual)  +  strata(ID_strata),
>method = "efron", data = new_data, x=T, y=T)
>  coefficients(mod)
>}
> 
>boot_data <- boot(data=data, statistic=boot.clogit, R=5000)
> 
> I have attached an overview of my data set.

You probably tried to attach something but you failed to note the section in 
the listinfo or posting guide where the list owners describe the rules for 
attachments. I think you would need to describe the sampling design more 
thoroughly. A simple description of the data layout may not be sufficient.

 The fact that you are clustering on individuals suggests you have some sort of 
repeated measures design and that you have somehow matched the individual to 
controls in some unstated ratio (handled by the strata. (Admittedly all 
guesswork and the more knowledgeable respondents (among which I'm not likely to 
reside)  are often hesitant to contribute substantive commentary unless they 
can narrow down range of possible design issues. I read Davison and Hinkley as 
suggesting that sampling by group but then keeping sampled groups undisturbed 
may have better chance of resulting in estimates of variances that match the 
superpopulation. See pages 100-102 of their book.

If my reading of that section is correct then I should think you would arrange 
you data so groups are in the long direction and single groups occupy a line of 
data with a single index. Then you would probably rearrange the data within the 
boot.clogit function so that the "inner" clogit call can handle it correctly.

-- 
David.
> 
> Thank you very much for your time.
> Best regards,
> Nell
> 
> 
> 
> 
> 
> 
>   [[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.

David Winsemius
Alameda, CA, USA

'Any technology distinguishable from magic is insufficiently advanced.'   
-Gehm's Corollary to Clarke's Third Law

__
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] error message for function: lmer (from lme4 package)

[David Winsemius]

> On Nov 14, 2017, at 5:13 AM, Fix Ace via R-help  wrote:
> 
> Dear R Community,
> My data have 3 conditions and each condition has 6 replicates. I am trying to 
> fit my data for a linear mixed model using the lmer function from lme4 
> package to find the random effects of the replicates;

Better venue for this question might be SIG-mixed-models. See the link 
avaialble at the bottom of every posting from rhelp:

https://stat.ethz.ch/mailman/listinfo/r-help:



> however, I got the error message. Here are the example codes:
>> example.3=data.frame(levels=as.numeric(XXX[,c(4)]),replicate=rep(c("0","1","2","3","4","5"),3),conditions=c(rep("11",6),rep("12",6),rep("13",6)))>
>>  example.3levels replicate conditions1  43. 0 112  
>> 42.0942 1 113  57.8131 2 114  57.1726
>>  3 115  77.8678 4 116  44.7578 5 117 
>>  69.5078 0 128  52.0581 1 129  40.0602   
>>   2 1210 45.5487 3 1211 43.6201 4 
>> 1212 60.4939 5 1213 64.1932 0 1314 53.4055   
>>   1 1315 59.6701 2 1316 52.6922 3
>>  1317 53.8712 4 1318 60.2770 5 13> 
>> m.example.3=lmer(as.numeric(levels)~conditions+(conditions|replicate),data=example.3)Error:
>>  number of observations (=18) <= number of random effects (=18) for term 
>> (conditions | replicate); the random-effects parameters and the residual 
>> variance (o
 r scale parameter) are probably unidentifiable> 
> Could anyone help me figure out how to fix the issue? 
> Thank you very much for any inputs!
> Ace
> 
>   [[alternative HTML version deleted]]

Complete mess. If you haven't yet been advised to posting in plain text, then 
this should be your wakeup call. If you have, then why are you ignoring 
sensible advice?


> 
> __
> 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.

David Winsemius
Alameda, CA, USA

'Any technology distinguishable from magic is insufficiently advanced.'   
-Gehm's Corollary to Clarke's Third Law

__
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] error message for function: lmer (from lme4 package)

[David Winsemius]

> On Nov 14, 2017, at 12:49 PM, Fix Ace  wrote:
> 
> Hi, David,
> 
> Thank you very much for getting back to me! Sorry about the messy code 
> example. I am re-posting here (including the error message):
> 
> > example.3=data.frame(levels=as.numeric(XXX[,c(4)]),replicate=rep(c("0","1","2","3","4","5"),3),conditions=c(rep("11",6),rep("12",6),rep("13",6)))
> > example.3
> levels replicate conditions
> 1  43. 0 11
> 2  42.0942 1 11
> 3  57.8131 2 11
> 4  57.1726 3 11
> 5  77.8678 4 11
> 6  44.7578 5 11
> 7  69.5078 0 12
> 8  52.0581 1 12
> 9  40.0602 2 12
> 10 45.5487 3 12
> 11 43.6201 4 12
> 12 60.4939 5 12
> 13 64.1932 0 13
> 14 53.4055 1 13
> 15 59.6701 2 13
> 16 52.6922 3 13
> 17 53.8712 4 13
> 18 60.2770 5 13
> > m.example.3=lmer(as.numeric(levels)~conditions+(conditions|replicate),data=example.3)
> Error: number of observations (=18) <= number of random effects (=18) for 
> term (conditions | replicate); the random-effects parameters and the residual 
> variance (or scale parameter) are probably unidentifiable

The error message seems fairly clear. The formula you have provided is asking 
for estimation of too many parameters. I think you probably want:

m.example.3=lmer(levels~conditions+(1|replicate),data=example.3)

... although your description of the hypothesis under test is ... non-existent.

-- 
David.
> > 
> 
> Please let me know if it is readable this time. 
> 
> Again, many thanks for your time and please help me fix the issue.
> 
> Kind regards,
> 
> Ace
> 
> 
> On Tuesday, November 14, 2017 12:19 PM, David Winsemius 
>  wrote:
> 
> 
> 
> > On Nov 14, 2017, at 5:13 AM, Fix Ace via R-help  
> > wrote:
> > 
> > Dear R Community,
> > My data have 3 conditions and each condition has 6 replicates. I am trying 
> > to fit my data for a linear mixed model using the lmer function from lme4 
> > package to find the random effects of the replicates;
> 
> Better venue for this question might be SIG-mixed-models. See the link 
> avaialble at the bottom of every posting from rhelp:
> 
> https://stat.ethz.ch/mailman/listinfo/r-help:
> 
> 
> 
> 
> > however, I got the error message. Here are the example codes:
> >> example.3=data.frame(levels=as.numeric(XXX[,c(4)]),replicate=rep(c("0","1","2","3","4","5"),3),conditions=c(rep("11",6),rep("12",6),rep("13",6)))>
> >>  example.3levels replicate conditions1  43.0 112  
> >> 42.09421113  57.81312114  57.17263 
> >>115  77.86784116  44.75785117  
> >> 69.50780128  52.05811129  40.06022 
> >>1210 45.548731211 43.620141212 
> >> 60.493951213 64.193201314 53.40551 
> >>1315 59.670121316 52.692231317 
> >> 53.871241318 60.2770513> 
> >> m.example.3=lmer(as.numeric(levels)~conditions+(conditions|replicate),data=example.3)Error:
> >>  number of observations (=18) <= number of random effects (=18) for term 
> >> (conditions | replicate); the random-effects parameters and the residual 
> >> variance (or scale parameter) are probably u
 nidentifiable> 
> > Could anyone help me figure out how to fix the issue? 
> > Thank you very much for any inputs!
> > Ace
> 
> > 
> > [[alternative HTML version deleted]]
> 
> Complete mess. If you haven't yet been advised to posting in plain text, then 
> this should be your wakeup call. If you have, then why are you ignoring 
> sensible advice?
> 
> 
> > 
> > __
> > 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.
> 
> David Winsemius
> Alameda, CA, USA
> 
> 'Any technology distinguishable from magic is insufficiently advanced.'  
> -Gehm's Corollary to Clarke's Third Law
> 
> 
> 
> 
> 
> 
> 

David Winsemius
Alameda, CA, USA

'Any technology distinguishable from magic is insufficiently advanced.'   
-Gehm's Corollary to Clarke's Third Law

__
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] Problems installing mice package

[David Winsemius]

> On Nov 15, 2017, at 1:08 AM, Jeremie Juste  wrote:
> 
> 
> 
> Hello,
> 
> I tried intalling mice package and got the following error:
> 
> * installing *source* package ‘mice’ ...
> ** package ‘mice’ successfully unpacked and MD5 sums checked
> ** libs
> g++  -I/usr/local/lib/R/include -DNDEBUG  
> -I"/home/djj/R/x86_64-pc-linux-gnu-library/3.4/Rcpp/include" 
> -I/usr/local/include   -fpic  -g -O2  -c RcppExports.cpp -o RcppExports.o
> g++  -I/usr/local/lib/R/include -DNDEBUG  
> -I"/home/djj/R/x86_64-pc-linux-gnu-library/3.4/Rcpp/include" 
> -I/usr/local/include   -fpic  -g -O2  -c match.cpp -o match.o
> g++ -shared -L/usr/local/lib -o mice.so RcppExports.o match.o Welcome to R! 
> Goodbye!

Somehow you have sent an R "welcome message" to the installer script. What was 
the action that started all this? Were you in an R session or was this 
something done from a bash console?

> g++: error: Welcome: No such file or directory
> g++: error: to: No such file or directory
> g++: error: R!: No such file or directory
> g++: error: Goodbye!: No such file or directory
> /usr/local/lib/R/share/make/shlib.mk:6: recipe for target 'mice.so' failed
> make: *** [mice.so] Error 1
> ERROR: compilation failed for package ‘mice’
> * removing ‘/home/djj/R/x86_64-pc-linux-gnu-library/3.4/mice’
> 
> The downloaded source packages are in
>   ‘/tmp/Rtmpgam70t/downloaded_packages’
> Error in library(mice) : there is no package called ‘mice’
> In addition: Warning message:
> In  :

That's a rather strange method for invoking install.packages. Usually it would 
be just:

install.packages("mice", repos = "http://cran.us.r-project.org";)



Your best venue for linux installation questions is probably:

https://stat.ethz.ch/mailman/listinfo/r-sig-db


>  installation of package ‘mice’ had non-zero exit status
> 
> I'm unable to resolve it. Any help please?
> 
> Best regards,
> 
> Jeremie
> 
> __
> 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.

David Winsemius
Alameda, CA, USA

'Any technology distinguishable from magic is insufficiently advanced.'   
-Gehm's Corollary to Clarke's Third Law

__
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] Autologistic regression in R

[David Winsemius]

> On Nov 14, 2017, at 4:39 PM, Mingke Li  wrote:
> 
> Hi,
> 
> I am new to autologistic regression and R. I do have questions when starting 
> a project in which I believe autologistic regression is needed.
> 
> I have a point layer whose attribute table stores the values of the dependent 
> variable and all the independent variables. I hope to to fit an autologistic 
> model to analyze which factors or combinations of factors have effects on the 
> presence/absence of the dependent variable (1 or 0).
> 
> I found other papers which applied autologistic regression in their study 
> almost used a grid system and defined their window sizes. So, my question is 
> do I have to convert my point layer to a grid system if I want to do this 
> analysis with R?
> 
> Also, what should I consider when I generate the grid system? How to 
> determine a proper size of cells? How about the searching window sizes?

Have you read the Posting Guide?

-- 
David.
> 
> Many Thanks.
> 
> Erin
> 
> 
>   [[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.

David Winsemius
Alameda, CA, USA

'Any technology distinguishable from magic is insufficiently advanced.'   
-Gehm's Corollary to Clarke's Third Law

__
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] error message for function: lmer (from lme4 package)

[David Winsemius]

> snipped
> 
> Hi, Bert,
> Sorry about that! David seemed to be able to read the post since he replied.

The only reason it appear readable was that you copied me as well as the list. 
Then HTML was still  the basic format although it did npt appear that way to me 
since my reader could handle it without difficulty. Everybody else saw a 
"complete mess" as Bert accurately put it.

Again, learn to post in plain text. Learn to post data objects with dput.

-- 
David.
> Here I just email you the sample code and error message:
>> example.3=data.frame(levels= as.numeric(XXX[,c(4)]), 
>> replicate=rep(c("0","1","2"," 3","4","5"),3),conditions=c( 
>> rep("11",6),rep("12",6),rep(> example.3levels replicate conditions1  
>> 43. 0 112  42.0942 1 113  57.8131
>>  2 114  57.1726 3 115  77.8678 4 116 
>>  44.7578 5 117  69.5078 0 128  52.0581   
>>   1 129  40.0602 2 1210 45.5487 3 
>> 1211 43.6201 4 1212 60.4939 5 1213 64.1932   
>>   0 1314 53.4055 1 1315 59.6701 2
>>  1316 52.6922 3 1317 53.8712 4 1318 60.2770  
>>5 13> m.example.3=lmer(as.numeric( levels)~conditions+( 
>> conditions|replicate),data= example.3)Error in lmer(as.numeric(levels) ~ 
>> conditions + (conditions | replicate),  :   could not find function "lmer"> 
> Hopefully you could read it and provide some comments!
> Thank you very much for your time.
> Kind regards,
> Ace

snipped

> 
>   [[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.

David Winsemius
Alameda, CA, USA

'Any technology distinguishable from magic is insufficiently advanced.'   
-Gehm's Corollary to Clarke's Third Law

__
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] family

[David Winsemius]

> On Nov 17, 2017, at 4:28 PM, Val  wrote:
> 
> Hi all,
> I am reading a huge data set(12M rows) that contains family information,
> Offspring, Parent1 and Parent2
> 
> Parent1 and parent2 should be in the first column as an offspring
> before their offspring information. Their parent information (parent1
> and parent2) should be  set to zero, if unknown.  Also the first
> column should be unique.
> 
> 
> Here is my sample data  set  and desired output.
> 
> 
> fam <- read.table(textConnection(" offspring  Parent1 Parent2
> Smith Alex1  Alexa
> Carla Alex1 0
> Jacky Smith   Abbot
> Jack  0   Jacky
> Almo  JackCarla
> "),header = TRUE)
> 
> 
> 
> desired output.
> Offspring Parent1 Parent2
> Alex1  00
> Alexa  00
> Abbot  00
> SmithAlex1  Alexa
> CarlaAlex1  0
> JackySmith   Abbot
> Jack   0 Jacky
> Almo JackCarla

You might get useful ideas by looking at ?'%in%" and ?union (set operations)

> fam$Parent1[!fam$Parent1 %in% fam$offspring]
[1] "Alex1" "Alex1" "0"
> fam$Parent2[!fam$Parent1 %in% fam$offspring]
[1] "Alexa" "0" "Jacky"

David.
> 
> Thank you.
> 
> __
> 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.

David Winsemius
Alameda, CA, USA

'Any technology distinguishable from magic is insufficiently advanced.'   
-Gehm's Corollary to Clarke's Third Law

__
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] Complicated analysis for huge databases

[David Winsemius]

s 
>>> combinations, this is followed by another group of customers with similar 
>>> meals combinations but different from the first group and so on. The 
>>> dataset looks like this :-
>>> 
>>> 
>>>> MyData
>>> 
>>>  Meal A Meal B Cust.ID  IIIIII IV   
>>> .. 600
>>> 
>>> 133 55 1 0   12 
>>>   0
>>> 
>>> 233 55  3 1  02 
>>>2
>>> 
>>> 333 55  5 2  11 
>>> 2
>>> 
>>> 444 66   70  2 
>>> 22
>>> 
>>> 5   44  66   41  1  
>>> 0   1
>>> 
>>> 6   44  6692  0 
>>>  1   2
>>> 
>>> .
>>> 
>>> .
>>> 
>>> 600,000
>>> 
>>> 
>>> 
>>> I wanted to find maf() for each column(from 4 to 600) after calculating the 
>>> frequency of the 3 values (0,1,2) but this should be done group by group 
>>> (i.e. group(33-55) : rows 1:3 then group(44-66) :rows 4:6 and so on).
>>> 
>>> 
>>> I can do the analysis  for the entire column but not group by group like 
>>> this :
>>> 
>>> 
>>> MAF <- apply(MyData[,4:600], 2, function(x)maf(tabulate(x+1)))
>>> 
>>> How can I modify this code to tell R to do the analysis group by group for 
>>> each column so I get maf value for 33-55 group of clolumn I, then maf value 
>>> for group 44-66 in the same column I,then the rest of groups in this column 
>>> and do the same for the remaining columns.
>>> 
>>> In fact, I'm interested in doing this analysis for only 300 columns but all 
>>> of the 600 columns.
>>> I have another sheet contains names of columns of interest like this :
>>> 
>>>> ColOfinterest
>>> 
>>> Col
>>> I
>>> IV
>>> V
>>> .
>>> .
>>> 300
>>> 
>>> Any one would help with the best combination of syntax to perform this 
>>> complex analysis?
>>> 
>>> Regards
>>> Allaisone
>>> 
>>> 
>>> 
>>> 
>>> 
>>> 
>>> 
>>>   [[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.

David Winsemius
Alameda, CA, USA

'Any technology distinguishable from magic is insufficiently advanced.'   
-Gehm's Corollary to Clarke's Third Law

__
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] help

[David Winsemius]

> On Nov 21, 2017, at 8:15 AM, yadav neog  wrote:
> 
> thank you for your valuable reply. I  have attached my commands, results, and
> data with this mail..maybe it will be beneficial for you to feedback.

If you had read the Posting Guide (and lsitinfo), you should have noted the 
special requirements regarding  attachments. The only person who received your 
attachments was Jeff. So you are essentially at his mercy regarding how he 
handles private communications.

-- 
David.


> 
> On Tue, Nov 21, 2017 at 9:13 PM, Jeff Newmiller 
> wrote:
> 
>> Your example is incomplete... as the bottom of this and every post says,
>> we need to be able to proceed from an empty R environment to wherever you
>> are having the problem (reproducible), in as few steps as possible
>> (minimal). The example needs to include data, preferably in R syntax as the
>> dput function creates... see the howtos referenced below for help with
>> that.  [1], [2], [3]
>> 
>> You also need to set your email program to send plain text format, since
>> HTML gets mangled to various degrees as it gets forced into text format
>> going through the mailing list. Read the Posting Guide.
>> 
>> A wild guess is that you have negative values in your data or too few data
>> points...
>> 
>> [1] http://stackoverflow.com/questions/5963269/how-to-make-
>> a-great-r-reproducible-example
>> 
>> [2] http://adv-r.had.co.nz/Reproducibility.html
>> 
>> [3] https://cran.r-project.org/web/packages/reprex/index.html (read the
>> vignette)
>> --
>> Sent from my phone. Please excuse my brevity.
>> 
>> On November 21, 2017 12:48:08 AM PST, yadav neog 
>> wrote:
>>> I am working on  Johansen cointegration test,  using urca and var
>>> package.
>>> in the selection of var, I have got following results.
>>> 
>>>> VARselect(newd, lag.max = 10,type = "none")
>>> 
>>> $selection
>>> AIC(n)  HQ(n)  SC(n) FPE(n)
>>>6  6  6  5
>>> 
>>> $criteria
>>>  1 2 3 4
>>> 56789
>>> AIC(n) -3.818646e+01 -3.864064e+01 -3.833435e+01 -4.089169e+01
>>> NaN -Inf -Inf -Inf -Inf
>>> HQ(n)  -3.754345e+01 -3.744647e+01 -3.658903e+01 -3.859523e+01
>>> NaN -Inf -Inf -Inf -Inf
>>> SC(n)  -3.630096e+01 -3.513899e+01 -3.321655e+01 -3.415775e+01
>>> NaN -Inf -Inf -Inf -Inf
>>> FPE(n)  2.700145e-17  2.114513e-17  5.350381e-17  2.035215e-17
>>> -9.147714e-650000
>>>10
>>> AIC(n) -Inf
>>> HQ(n)  -Inf
>>> SC(n)  -Inf
>>> FPE(n)0
>>> 
>>> Warning messages:
>>> 1: In log(sigma.det) : NaNs produced
>>> 2: In log(sigma.det) : NaNs produced
>>> 3: In log(sigma.det) : NaNs produced
>>> 
>>> so do in my ca.jo test... I have found similar NaNs results. please
>>> help me
>>> in solving the problem.
>>> Yadawananda Neog
>>> Research Scholar
>>> Department of Economics
>>> Banaras Hindu University
>>> Mob. 9838545073
>>> 
>>>  [[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.
>> 
> 
> 
> 
> -- 
> Yadawananda Neog
> Research Scholar
> Department of Economics
> Banaras Hindu University
> Mob. 9838545073
> __
> 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.

David Winsemius
Alameda, CA, USA

'Any technology distinguishable from magic is insufficiently advanced.'   
-Gehm's Corollary to Clarke's Third Law

__
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] R-How to unlist data frame multiple structured list column value and new column

[David Winsemius]

quot;st7, st9"), 
> .Names = c("PassID",
>"Name", "Name1", "PassIt"), class = "data.frame", row.names = 1L)),
>TrainnerDe = list(structure(list(Trnid = 1, Trncont = 5, EmpAddInfo = list(
>structure(list(CohID = "pi", InVoice = 77), .Names = c("CohID",
>"InVoice"), class = "data.frame", row.names = 1L))), .Names = 
> c("Trnid",
>"Trncont", "EmpAddInfo"), class = "data.frame", row.names = 1L),
>structure(list(Trnid =7,Trncont = 3, EmpAddInfo = list(structure(list(
>CohID = c("AB", "NI", "OL"), InVoice = c("4", "Y")), .Names = 
> c("CohID",
>"InVoice"), class = "data.frame", row.names = 1L))), .Names = 
> c("Trnid",
>"Trncont", "EmpAddInfo"), class = "data.frame", row.names = 1L),
> structure(list(Trnid =c("U","l"),Trncont =c("10","78"), EmpAddInfo = 
> list(structure(list(
>CohID = c("AB", "NI", "OL"), InVoice = c("4", "Y")), .Names = 
> c("CohID",
>"InVoice"), class = "data.frame", row.names = 1L))), .Names = 
> c("Trnid",
>"Trncont", "EmpAddInfo"), class = "data.frame", row.names = 1L),
>structure(list(Trnid = c("1","3"), Trncont = c("yt","re"), EmpAddInfo 
> = list(structure(list(
>CohID = c("Ly","qp"), InVoice = c("2","P")), .Names = c("CohID", 
> "InVoice"
>), class = "data.frame", row.names = 1L))), .Names = c("Trnid",
>"Trncont", "EmpAddInfo"), class = "data.frame", row.names = 1L))), 
> .Names = c("ID",
> "ContractDe", "PassengersDe", "TrainnerDe"), row.names = c(NA, 4L), class = 
> "data.frame")
> 
> 
> i was stuck on this change part of my project, please help me to out this. 
> Thanks.
> 
> 
> 
> Sent from Outlook<http://aka.ms/weboutlook>
> 
>   [[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.

David Winsemius
Alameda, CA, USA

'Any technology distinguishable from magic is insufficiently advanced.'   
-Gehm's Corollary to Clarke's Third Law

__
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] libPaths displays truncated path?

[David Winsemius]

> On Nov 23, 2017, at 4:34 AM, Loris Bennett  wrote:
> 
> Hi,
> 
> TL;DR
> -
> 
>  I define the path
> 
>/cm/shared/apps/R/site-library/3.4.2
> 
>  and add it to libPath.  Why does libPath then display it as
> 
>/cm/shared/apps/R/site-library/3.4

Generally one only has a different library for each major version of R. Major 
versions are consider just the first two numbers in the dot-separated versions 
system. Apparently libPaths is "smart" enough to make an effort to adhere to 
that convention.  It appears your definition of "major" differs from the usual 
convention.

-- 
David
> 
>  ?
> 
> Long version
> 
> 
> I run a cluster of diskless nodes for which the OS is loaded
> directly into RAM and other software is provided by an NFS server.
> However, in the case of R, we use the R version provided by the OS and
> just install additional packages on the NFS server.
> 
> So that R can find these additional packages, I have the following in
> the site-wide Rprofile
> 
>  v <- R.Version()
>  base_path = "/cm/shared/apps/R/site-library/"
>  major_minor_version = paste(v["major"],v["minor"],sep=".")
>  cm_shared_lib_path = paste0(base_path,major_minor_version)
>  full_cm_shared_lib_path <- c(file.path(chartr("\\", "/", R.home()), 
> "site-library"), cm_shared_lib_path)
>  .libPaths( c( .libPaths(), full_cm_shared_lib_path ) )
> 
> Thus, when I start R I get this:
> 
>> full_cm_shared_lib_path
>  [1] "/usr/lib64/R/site-library"   
>  [2] "/cm/shared/apps/R/site-library/3.4.2"
> 
> but also this
> 
>> .libPaths()
>  [1] "/home/loris/R/x86_64-redhat-linux-gnu-library/3.4"
>  [2] "/usr/lib64/R/library" 
>  [3] "/usr/share/R/library" 
>  [4] "/usr/lib64/R/site-library"
>  [5] "/cm/shared/apps/R/site-library/3.4"   
> 
> However, in order to get R to find the packages, I have to add a
> symbolic link, thus:
> 
>  [/cm/shared/apps/R/site-library] $ ls -l 3.4.2
>  lrwxrwxrwx 1 root root 3 Nov 23 09:21 3.4.2 -> 3.4
> 
> So, my mistake was to think that "minor" would return "4", whereas it in
> fact returns "4.2".  So I actually set the path to ".../3.4.2" and
> that's where R looks for packages.
> 
> But why does libPaths display the the path I *thought* I had set, but in
> fact looks at the path I *really did* set?
> 
> Cheers,
> 
> Loris
> 
> -- 
> Dr. Loris Bennett (Mr.)
> ZEDAT, Freie Universität Berlin Email loris.benn...@fu-berlin.de
> 
> __
> 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.

David Winsemius
Alameda, CA, USA

'Any technology distinguishable from magic is insufficiently advanced.'   
-Gehm's Corollary to Clarke's Third Law

__
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 sets of proportion tests

[David Winsemius]

ns and have only the columns with 
>> significant p-value results to be written in the the third row under each 
>> column so the final output has to be something like this :-
>> 
>> 
>>  variable IVariable III  
>> .
>> 
>> Freq.of cases (gp I)  6493   5524
>> 
>> Freq. of cases (gpII) 509  54
>> 
>> p-values  0.02   0.010
>> 
>> Note, for example, that the 2nd column has bee removed as it resulted in a 
>> non-significant p-value result while col 1 and col 3 were included since 
>> p-value is less than 0.05.
>> 
>> I'm not sure how to get the p-values only without other details but for the 
>> analysis itself , I believe it can be done with apply() function but its not 
>> clear to me how to specify the 2nd argument(n=samlpe sizes) in the prop.test.
>> 
>> MyResults <- apply(Mydata, 2, function(x)prop.test(Mydata,c(200,100))
>> 
>> How can I modify the "n" argument part to solve the issue of non-equivalent 
>> length between "x" and "n" ?. How can I modify this further to return only 
>> significant p-values results ?. Any help would be very appreciated ..
>> 
>> Regards
>> 
>>[[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.

David Winsemius
Alameda, CA, USA

'Any technology distinguishable from magic is insufficiently advanced.'   
-Gehm's Corollary to Clarke's Third Law

__
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] dplyr - add/expand rows

[David Winsemius]

> On Nov 25, 2017, at 11:18 AM, Hutchinson, David (EC) 
>  wrote:
> 
> I have a returned tibble of station operational record similar to the 
> following:
> 
>> data.collection
> # A tibble: 5 x 4
>  STATION_NUMBER YEAR_FROM YEAR_TO RECORD
> 
> 107EA001  19601960QMS
> 207EA001  19611970QMC
> 307EA001  19711971QMM
> 407EA001  19721976QMC
> 507EA001  19771983QRC
> 
> I would like to reshape this to one operational record (row) per year per 
> station. Something like:
> 
> 07EA001  1960  QMS
> 07EA001  1961  QMC
> 07EA001  1962  QMC
> 07EA001  1963  QMC
> ...
> 07EA001  1971  QMM
> 
> Can this be done in dplyr easily?

Probably, yes. This looks like a feasible plan might be to "fill-in" the gaps 
with a last observation carried forward value within categories of station 
number. The na.locf function in package zoo is very handy for some of these 
tasks. Or Perhaps merging this data with a skeleton data object with 
station numbers and a `seq`-built vectors for the range of years.  Why don't 
you post a data example with sufficient complexity to represent the problem? 
Perhaps:

 dput( head( data.collection, 20) )

It's clear that the first 5 lines are not sufficient since there's only one 
station. It's kind of a pain to try to construct tibble objects from their 
print output representations. And posting code to build examples is a specific 
suggestion in the Posting Guide.

-- 
David
> 
> Thanks in advance,
> 
> David
> 
>   [[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.

David Winsemius
Alameda, CA, USA

'Any technology distinguishable from magic is insufficiently advanced.'   
-Gehm's Corollary to Clarke's Third Law

__
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] 2^3 confounded factorial experiment

[David Winsemius]

> On Nov 29, 2017, at 9:20 AM, Jyoti Bhogal  wrote:
> 
> The following R commands were written:
>> help.search("factorial")
>> data(npk)
>> npk
>> coef(npk.aov)
> 
> In the output of coef command, please explain me the interpretation of 
> coefficients of block1 to block 6 in this 2^3 confounded factorial experiment.

This is very much a statistics question and as such is off-topic (as it also 
would be off-topic on StackOverflow.) Rhelp is for persons having difficulty 
coding the R language itself.

Consider CrossValidated.com but read their posting help section first since 
this is really very terse question. Better would be to include the output and 
make your best interpretation so peolple get the sense you at least put in some 
individual effort.

> 
> Thanks.
>   [[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.

David Winsemius
Alameda, CA, USA

'Any technology distinguishable from magic is insufficiently advanced.'   
-Gehm's Corollary to Clarke's Third Law

__
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] Help avoiding setting column type two times

[David Winsemius]

> On Nov 30, 2017, at 6:27 AM, Martin Møller Skarbiniks Pedersen 
>  wrote:
> 
> Hi,
>  I think and hope this a good place to ask for code review for a R
> beginners?
> 
>  I have made a R script which generates a dataset based on 2009 danish
> referendum and it does work.
> 
>  But I think the code could be better and I would any comments how the
> code can be improved.
>  At least I would like to know how I avoid converting several of the
> columns to factors in the end of the code?
> 
> Description of the code:
> 
>  It reads a lot of xml-files from ../raw/ and saves a data.frame with
> information
> from these xml-files.
> 
>  In the ../raw/ directiory I have placed the xml-files which I got from
> "Statistics Denmark"
>  I have also put these xml-files on my website and they can be download
> freely from http://20dage.dk/R/referendum-2009/raw.tar.gz
> 
>  The code is below but I have also put the code at this place:
> http://20dage.dk/R/referendum-2009/convert_from_xml.R
> 
> Best Regards
> Martin M. S. Pedersen
> 
> ---
> library(xml2)
> 
> convert_one_file <- function(url) {
>x <- read_xml(url)
> 
>Sted <- xml_find_first(x, ".//Sted")
>StedType <- xml_attr(Sted, "Type")
>StedTekst <- xml_text(Sted)
> 
>Parti <- xml_find_all(x, ".//Parti")
>PartiId <- xml_attr(Parti, "Id")
>PartiBogstav <- xml_attr(Parti, "Bogstav")
>PartiNavn <- xml_attr(Parti, "Navn")
> 
> 
>StemmerAntal <- xml_attr(Parti, "StemmerAntal")
>Stemmeberettigede <- xml_integer(xml_find_first(x,
> ".//Stemmeberettigede"))
>DeltagelsePct <- xml_double(xml_find_first(x, ".//DeltagelsePct"))
>IAltGyldigeStemmer <- xml_integer(xml_find_first(x,
> ".//IAltGyldigeStemmer"))
>BlankeStemmer <- xml_integer(xml_find_first(x, ".//BlankeStemmer"))
>AndreUgyldigeStemmer <- xml_integer(xml_find_first(x,
> ".//AndreUgyldigeStemmer"))
> 
>data.frame(cbind(StedType, StedTekst, PartiId, PartiBogstav, PartiNavn,
> StemmerAntal, Stemmeberettigede, DeltagelsePct,
> IAltGyldigeStemmer,
>   BlankeStemmer, AndreUgyldigeStemmer), stringsAsFactors = FALSE)

The construction `data.frame(cbind( ...` is a serious source of potential 
error. The cbind coerces to matrix class which also then coerces to a single 
atomic class, either numeric or character. Factors loose all their meaning. 
Dates get messed up. Error ensues. Better would be:

 data.frame( StedType, StedTekst, PartiId, PartiBogstav, PartiNavn,
 StemmerAntal, Stemmeberettigede, DeltagelsePct,
 IAltGyldigeStemmer, BlankeStemmer, AndreUgyldigeStemmer,
 stringsAsFactors = FALSE)

-- 
David.


> }
> 
> raw_path <- "../raw"
> filenames <- dir(path = raw_path, pattern = "fintal_.*", full.names = T)
> 
> result <- data.frame(StedType = factor(),
> StedTekst = character(),
> PartiId   = factor(),
> PartiBogstav = factor(),
> PartiNavn= factor(),
> StemmerAntal = integer(),
> Stemmeberettigede = integer(),
> DeltagelsePct = numeric(),
> IAltGyldigeStemmer = integer(),
> BlankeStemmer = integer(),
> AndreUgyldigeStemmer = integer(),
> stringsAsFactors = FALSE)
> 
> for (i in 1:length(filenames)) {
>#cat(paste0(filenames[i],"\n"))
>returnCode <-  tryCatch({
>   result <- rbind(result, convert_one_file(filenames[i]))
>}, error = function(e) {
>   cat(paste0(filenames[i]," failed:\n",e,"\n"))
>})
> }
> 
> result$StedType <- as.factor(result$StedType)
> result$PartiId <- as.factor(result$PartiId)
> result$PartiBogstav <- as.factor(result$PartiBogstav)
> result$PartiNavn <- as.factor(result$PartiNavn)
> result$StemmerAntal <- as.integer(result$StemmerAntal)
> result$Stemmeberettigede <- as.integer(result$Stemmeberettigede)
> result$DeltagelsePct <- as.numeric(result$DeltagelsePct)
> result$IAltGyldigeStemmer <- as.integer(result$IAltGyldigeStemmer)
> result$BlankeStemmer <- as.integer(result$BlankeStemmer)
> result$AndreUgyldigeStemmer <- as.integer(result$AndreUgyldigeStemmer)
> str(result)
> save(result, file = "folkeafstemning2009.Rdata")
> 
>   [[alternative HTML version deleted]]
> 
> __
> R-help@r-project.org mailing 

Re: [R] Odd dates generated in Forecasts

[David Winsemius]

 "08/01/09", "08/01/15", "08/02/03", "08/02/08", "08/02/14",
> "08/03/13", "08/04/07", "08/04/12", "08/05/06", "08/06/05", "08/06/11",
> "08/06/16", "08/07/04", "08/07/10", "08/08/09", "08/08/15", "08/09/03",
> "08/09/08", "08/09/14", "08/10/13", "08/11/07", "08/11/12", "08/12/06",
> "08/13/05", "08/13/11", "08/13/16", "08/14/04", "08/14/10", "08/15/09",
> "08/15/15", "08/16/03", "08/16/08", "08/16/14", "08/17/13", "08/18/07",
> "08/18/12", "08/19/06", "08/20/05", "08/20/11", "08/20/16", "08/21/04",
> "08/21/10", "08/22/09", "08/22/15", "08/23/03", "08/23/08", "08/23/14",
> "08/24/13", "08/25/07", "08/25/12", "08/26/06", "08/27/05", "08/27/11",
> "08/27/16", "08/28/04", "08/28/10", "08/29/09", "08/29/15", "08/30/03",
> "08/30/08", "08/30/14", "08/31/13", "09/01/07", "09/01/12", "09/02/06",
> "09/03/05", "09/03/11", "09/03/16", "09/04/04", "09/04/10", "09/05/09",
> "09/05/15", "09/06/03", "09/06/08", "09/06/14", "09/07/13", "09/08/07",
> "09/08/12", "09/09/06", "09/10/05", "09/10/11", "09/10/16", "09/11/04",
> "09/11/10", "09/12/09", "09/12/15", "09/13/03", "09/13/08", "09/13/14",
> "09/14/13", "09/15/07", "09/15/12", "09/16/06", "09/17/05", "09/17/11",
> "09/17/16", "09/18/04", "09/18/10", "09/19/09", "09/19/15", "09/20/03",
> "09/20/08", "09/20/14", "09/21/13", "09/22/07", "09/22/12", "09/23/06",
> "09/24/05", "09/24/11", "09/24/16", "09/25/04", "09/25/10", "09/26/09",
> "09/26/15", "09/27/03", "09/27/08", "09/27/14", "09/28/13", "09/29/07",
> "09/29/12", "09/30/06", "10/01/05", "10/01/11", "10/01/16", "10/02/04",
> "10/02/10", "10/03/09", "10/03/15", "10/04/03", "10/04/08", "10/04/14",
> "10/05/13", "10/06/07", "10/06/12", "10/07/06", "10/08/05", "10/08/11",
> "10/08/16", "10/09/04", "10/09/10", "10/10/09", "10/10/15", "10/11/03",
> "10/11/08", "10/11/14", "10/12/13", "10/13/07", "10/13/12", "10/14/06",
> "10/15/05", "10/15/11", "10/15/16", "10/16/04", "10/16/10", "10/17/09",
> "10/17/15", "10/18/03", "10/18/08", "10/18/14", "10/19/13", "10/20/07",
> "10/20/12", "10/21/06", "10/22/05", "10/22/11", "10/22/16", "10/23/04",
> "10/23/10", "10/24/09", "10/24/15", "10/25/03", "10/25/08", "10/25/14",
> "10/26/13", "10/27/07", "10/27/12", "10/28/06", "10/29/05", "10/29/11",
> "10/29/16", "10/30/04", "10/30/10", "10/31/09", "10/31/15", "11/01/03",
> "11/01/08", "11/01/14", "11/02/13", "11/03/07", "11/03/12", "11/04/06",
> "11/05/05", "11/05/11", "11/05/16", "11/06/04", "11/06/10", "11/07/09",
> "11/07/15", "11/08/03", "11/08/08", "11/08/14", "11/09/13", "11/10/07",
> "11/10/12", "11/11/06", "11/12/05", "11/12/11", "11/12/16", "11/13/04",
> "11/13/10", "11/14/09", "11/14/15", "11/15/03", "11/15/08", "11/15/14",
> "11/16/13", "11/17/07", "11/17/12", "11/18/06",

Re: [R] Odd dates generated in Forecasts

[David Winsemius]

> On Dec 6, 2017, at 11:09 AM, Paul Bernal  wrote:
> 
> Thank you very much David. As a matter of fact, I solved it by doing the 
> following:
> 
> MyTimeSeriesObj <- ts(MyData, freq=365.25/7, 
> start=decimal_date(mdy("01-04-2003")))
> 
> After doing that adjustment, my forecasts dates started from 2017 on.

Not clear what MyData consisted of. If it was the factor variable in the 
`dataset` object in your earlier communication, then you gave the `ts`  
function garbage. The fact that it was labeled to your liking would not create 
a palatable result.

-- 
David.
> 
> Cheers,
> 
> Paul
> 
> 2017-12-06 12:03 GMT-05:00 David Winsemius :
> 
> > On Dec 6, 2017, at 5:07 AM, Paul Bernal  wrote:
> >
> > Dear friends,
> >
> > I have a weekly time series which starts on Jan 4th, 2003 and ends on
> > december 31st, 2016.
> >
> > I set up my ts object as follows:
> >
> > MyTseries <- ts(mydataset, start=2003, end=2016, frequency=52)
> >
> > MyModel <- auto.arima(MyTseries, d=1, D=1)
> >
> > MyModelForecast <- forecast (MyModel, h=12)
> >
> > Since my last observation was on december 31st, 2016 I expected my forecast
> > date to start on 2017, but instead it returned odd dates.
> >
> > This is my dataset
> >
> > dput(dataset):
> >
> > structure(list(Date = structure(c(8L, 22L, 36L, 50L, 64L, 78L,
> > 92L, 106L, 120L, 134L, 148L, 162L, 176L, 190L, 204L, 218L, 232L,
> > 246L, 260L, 274L, 288L, 302L, 316L, 330L, 344L, 358L, 372L, 386L,
> > 400L, 414L, 428L, 442L, 456L, 470L, 484L, 498L, 512L, 526L, 540L,
> > 554L, 568L, 582L, 596L, 610L, 624L, 638L, 652L, 666L, 680L, 694L,
> > 708L, 722L, 5L, 19L, 33L, 47L, 61L, 75L, 89L, 103L, 117L, 130L,
> > 144L, 158L, 172L, 186L, 200L, 214L, 228L, 242L, 256L, 270L, 284L,
> > 298L, 312L, 326L, 340L, 354L, 368L, 382L, 396L, 410L, 424L, 438L,
> > 452L, 466L, 480L, 494L, 508L, 522L, 536L, 550L, 564L, 578L, 592L,
> > 606L, 620L, 634L, 648L, 662L, 676L, 690L, 704L, 718L, 1L, 15L,
> > 29L, 43L, 57L, 71L, 85L, 99L, 113L, 127L, 141L, 155L, 169L, 183L,
> > 197L, 211L, 225L, 239L, 253L, 267L, 281L, 295L, 309L, 323L, 337L,
> > 351L, 365L, 379L, 393L, 407L, 421L, 435L, 449L, 463L, 477L, 491L,
> > 505L, 519L, 533L, 547L, 561L, 575L, 589L, 603L, 617L, 631L, 645L,
> > 659L, 673L, 687L, 701L, 715L, 729L, 13L, 27L, 41L, 55L, 69L,
> > 83L, 97L, 111L, 126L, 140L, 154L, 168L, 182L, 196L, 210L, 224L,
> > 238L, 252L, 266L, 280L, 294L, 308L, 322L, 336L, 350L, 364L, 378L,
> > 392L, 406L, 420L, 434L, 448L, 462L, 476L, 490L, 504L, 518L, 532L,
> > 546L, 560L, 574L, 588L, 602L, 616L, 630L, 644L, 658L, 672L, 686L,
> > 700L, 714L, 728L, 12L, 26L, 40L, 54L, 68L, 82L, 96L, 110L, 124L,
> > 138L, 152L, 166L, 180L, 194L, 208L, 222L, 236L, 250L, 264L, 278L,
> > 292L, 306L, 320L, 334L, 348L, 362L, 376L, 390L, 404L, 418L, 432L,
> > 446L, 460L, 474L, 488L, 502L, 516L, 530L, 544L, 558L, 572L, 586L,
> > 600L, 614L, 628L, 642L, 656L, 670L, 684L, 698L, 712L, 726L, 10L,
> > 24L, 38L, 52L, 66L, 80L, 94L, 108L, 121L, 135L, 149L, 163L, 177L,
> > 191L, 205L, 219L, 233L, 247L, 261L, 275L, 289L, 303L, 317L, 331L,
> > 345L, 359L, 373L, 387L, 401L, 415L, 429L, 443L, 457L, 471L, 485L,
> > 499L, 513L, 527L, 541L, 555L, 569L, 583L, 597L, 611L, 625L, 639L,
> > 653L, 667L, 681L, 695L, 709L, 723L, 6L, 20L, 34L, 48L, 62L, 76L,
> > 90L, 104L, 118L, 132L, 146L, 160L, 174L, 188L, 202L, 216L, 230L,
> > 244L, 258L, 272L, 286L, 300L, 314L, 328L, 342L, 356L, 370L, 384L,
> > 398L, 412L, 426L, 440L, 454L, 468L, 482L, 496L, 510L, 524L, 538L,
> > 552L, 566L, 580L, 594L, 608L, 622L, 636L, 650L, 664L, 678L, 692L,
> > 706L, 720L, 3L, 17L, 31L, 45L, 59L, 73L, 87L, 101L, 115L, 131L,
> > 145L, 159L, 173L, 187L, 201L, 215L, 229L, 243L, 257L, 271L, 285L,
> > 299L, 313L, 327L, 341L, 355L, 369L, 383L, 397L, 411L, 425L, 439L,
> > 453L, 467L, 481L, 495L, 509L, 523L, 537L, 551L, 565L, 579L, 593L,
> > 607L, 621L, 635L, 649L, 663L, 677L, 691L, 705L, 719L, 2L, 16L,
> > 30L, 44L, 58L, 72L, 86L, 100L, 114L, 128L, 142L, 156L, 170L,
> > 184L, 198L, 212L, 226L, 240L, 254L, 268L, 282L, 296L, 310L, 324L,
> > 338L, 352L, 366L, 380L, 394L, 408L, 422L, 436L, 450L, 464L, 478L,
> > 492L, 506L, 520L, 534L, 548L, 562L, 576L, 590L, 604L, 618L, 632L,
> > 646L, 660L, 674L, 688L, 702L, 716L, 730L, 14L, 28L, 42L, 56L,
> > 70L, 84L, 98L, 112L, 125L, 139L, 153L, 167L, 181L, 195L, 209L,
> > 223L, 237L, 251L, 265L, 279L, 293L, 307L, 321L, 335L, 349L, 363L,
> > 377L, 391L, 405L, 419L, 433L, 447L, 461L, 475L, 489L, 503L, 517L,
> > 531L, 545L, 559L, 573L, 587L, 601L, 615L, 629L, 643L, 657L, 671L,
> > 685L, 6

Re: [R] Remove

[David Winsemius]

> On Dec 6, 2017, at 3:15 PM, Ashta  wrote:
> 
> Hi all,
> In a data set I have group(GR) and two variables   x and y. I want to
> remove a  group that have  the same record for the x variable in each
> row.
> 
> DM <- read.table( text='GR x y
> A 25 125
> A 23 135
> A 14 145
> A 12 230
> B 25 321
> B 25 512
> B 25 123
> B 25 451
> C 11 521
> C 14 235
> C 15 258
> C 10 654',header = TRUE, stringsAsFactors = FALSE)
> 
> In this example the output should contain group A and C  as group B
> has   the same record  for the variable x .
> 
> The result will be
> A 25 125
> A 23 135
> A 14 145
> A 12 230
> C 11 521
> C 14 235
> C 15 258
> C 10 654

Try:

DM[ !duplicated(DM$x) , ]
> 
> How do I do it R?
> Thank you.
> 
> __
> 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.

David Winsemius
Alameda, CA, USA

'Any technology distinguishable from magic is insufficiently advanced.'   
-Gehm's Corollary to Clarke's Third Law

__
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] Remove

[David Winsemius]

> On Dec 6, 2017, at 4:27 PM, Ashta  wrote:
> 
> Thank you Ista! Worked fine.

Here's another (possibly more direct in its logic?):

 DM[ !ave(DM$x, DM$GR, FUN= function(x) {!length(unique(x))==1}), ]
  GR  x   y
5  B 25 321
6  B 25 512
7  B 25 123
8  B 25 451

-- 
David

> On Wed, Dec 6, 2017 at 5:59 PM, Ista Zahn  wrote:
>> Hi Ashta,
>> 
>> There are many ways to do it. Here is one:
>> 
>> vars <- sapply(split(DM$x, DM$GR), var)
>> DM[DM$GR %in% names(vars[vars > 0]), ]
>> 
>> Best
>> Ista
>> 
>> On Wed, Dec 6, 2017 at 6:58 PM, Ashta  wrote:
>>> Thank you Jeff,
>>> 
>>> subset( DM, "B" != x ), this works if I know the group only.
>>> But if I don't know that group in this case "B", how do I identify
>>> group(s) that  all elements of x have the same value?
>>> 
>>> On Wed, Dec 6, 2017 at 5:48 PM, Jeff Newmiller  
>>> wrote:
>>>> subset( DM, "B" != x )
>>>> 
>>>> This is covered in the Introduction to R document that comes with R.
>>>> --
>>>> Sent from my phone. Please excuse my brevity.
>>>> 
>>>> On December 6, 2017 3:21:12 PM PST, David Winsemius 
>>>>  wrote:
>>>>> 
>>>>>> On Dec 6, 2017, at 3:15 PM, Ashta  wrote:
>>>>>> 
>>>>>> Hi all,
>>>>>> In a data set I have group(GR) and two variables   x and y. I want to
>>>>>> remove a  group that have  the same record for the x variable in each
>>>>>> row.
>>>>>> 
>>>>>> DM <- read.table( text='GR x y
>>>>>> A 25 125
>>>>>> A 23 135
>>>>>> A 14 145
>>>>>> A 12 230
>>>>>> B 25 321
>>>>>> B 25 512
>>>>>> B 25 123
>>>>>> B 25 451
>>>>>> C 11 521
>>>>>> C 14 235
>>>>>> C 15 258
>>>>>> C 10 654',header = TRUE, stringsAsFactors = FALSE)
>>>>>> 
>>>>>> In this example the output should contain group A and C  as group B
>>>>>> has   the same record  for the variable x .
>>>>>> 
>>>>>> The result will be
>>>>>> A 25 125
>>>>>> A 23 135
>>>>>> A 14 145
>>>>>> A 12 230
>>>>>> C 11 521
>>>>>> C 14 235
>>>>>> C 15 258
>>>>>> C 10 654
>>>>> 
>>>>> Try:
>>>>> 
>>>>> DM[ !duplicated(DM$x) , ]
>>>>>> 
>>>>>> How do I do it R?
>>>>>> Thank you.
>>>>>> 
>>>>>> __
>>>>>> 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.
>>>>> 
>>>>> David Winsemius
>>>>> Alameda, CA, USA
>>>>> 
>>>>> 'Any technology distinguishable from magic is insufficiently advanced.'
>>>>> -Gehm's Corollary to Clarke's Third Law
>>>>> 
>>>>> __
>>>>> 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.

David Winsemius
Alameda, CA, USA

'Any technology distinguishable from magic is insufficiently advanced.'   
-Gehm's Corollary to Clarke's Third Law

__
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] Error in loadNamespace

[David Winsemius]

> On Dec 7, 2017, at 10:47 AM, James Henson  wrote:
> 
> Hello R Community,
> 
> I inadvertently updated packages via R Studio when a package was open.  Now
> when R Studio is opened the message below appears in the console panel.
> 
> Error in loadNamespace(name) : there is no package called ‘yaml’
> 
> Error in loadNamespace(name) : there is no package called ‘yaml’
> 
> When running R code, so far the only function that has not worked is
> datatable () in the ‘DT’ package.  Removing the ‘DT’ package and
> reinstalling it had no effect.  A possible fix is to update R from an older
> version, which was not running when I made the mistake.  Is this a good
> idea?

There is a package on CRAN named "yaml". Have your tried installing it?
> 
> Thanks for your help.
> 
> James F. Henson
> 
>   [[alternative HTML version deleted]]

Rhelp is a plain text mailing list.

> 
> __
> 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.

David Winsemius
Alameda, CA, USA

'Any technology distinguishable from magic is insufficiently advanced.'   
-Gehm's Corollary to Clarke's Third Law

__
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] Seeking help with code

[David Winsemius]

> On Dec 7, 2017, at 2:51 PM, Stephanie Tsalwa  wrote:
> 
> Assuming the days of raining during half a year of all states(provinces) of a 
> country is normally distributed (mean=�, standard deviation=�) with sigma (�) 
> equals to 2. We now have 10 data points here: 26.64, 30.65, 31.27, 33.04, 
> 32.56, 29.10, 28.96, 26.44, 27.76, 32.27. Try to get the 95% level of CI for 
> �, using parametric Bootstrap method with bootstrap size B=8000.

This has the definite appearance of a homework assignment. It also has several 
other features which likewise suggest a failure to read the Posting Guide, a 
link to which is at the bottom of every posting sent out from the Rhelp 
mail-server.

You are advised to post credential/affiliation information and more context of 
the wider goals if my assumption of HW-status is incorrect.

-- 
David.

(And yes I do realize that I am also not posting any affiliation information, 
either, but I'm retired and am not asking for free advice, am I? Don't follow 
my bad example.)
> 
> my code - what am i doing wrong
> 
> #set sample size n, bootstrap size B
> n = 10
> b = 8000
> 
> 
> set a vector of days of rain into "drain"
> drain = c(26.64, 30.65, 31.27, 33.04, 32.56, 29.10, 28.96, 26.44, 27.76, 
> 32.27)
> 
> #calculate mean of the sample for days of rain
> mdr=mean(drain)
> mdr
> 
> #calculate the parameter of the exponential distribution
> lambdahat = 1.0/mdr
> lambdahat
> 
> #draw the bootstrap sample from Exponential
> x = rexp(n*b, lambdahat)
> x
> 
> bootstrapsample = matrix(x, nrow=n, ncol=b)
> bootstrapsample
> 
> # Compute the bootstrap lambdastar
> lambdastar = 1.0/colMeans(bootstrapsample)
> lambdastar
> 
> # Compute the differences
> deltastar = lambdastar - lambdahat
> deltastar
> 
> # Find the 0.05 and 0.95 quantile for deltastar
> d = quantile(deltastar, c(0.05,0.95))
> d
> 
> # Calculate the 95% confidence interval for lambda.
> ci = lambdahat - c(d[2], d[1])
> ci
> 
> 
>   [[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.

David Winsemius
Alameda, CA, USA

'Any technology distinguishable from magic is insufficiently advanced.'   
-Gehm's Corollary to Clarke's Third Law

__
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] Remove

[David Winsemius]

> On Dec 8, 2017, at 4:48 PM, Ashta  wrote:
> 
> Hi David, Ista and all,
> 
> I  have one related question  Within one group I want to keep records
> conditionally.
> example within
> group A I want keep rows that have  " x" values  ranged  between 15 and 30.
> group B I want keep rows that have  " x" values  ranged  between  40 and 50.
> group C I want keep rows that have  " x" values  ranged  between  60 and 75.

When you have a problem where there are multiple "parallel: parameters, the 
function to "reach for" is `mapply`. 

mapply( your_selection_func, group_vec, min_vec, max_vec)

... and this will probably return the values as a list (of dataframes if you 
build the function correctly,  so you may may need to then do:

do.call(rbind, ...)

-- 
David.
> 
> 
> DM <- read.table( text='GR x y
> A 25 125
> A 23 135
> A 14 145
> A 35 230
> B 45 321
> B 47 512
> B 53 123
> B 55 451
> C 61 521
> C 68 235
> C 85 258
> C 80 654',header = TRUE, stringsAsFactors = FALSE)
> 
> 
> The end result will be
> A 25 125
> A 23 135
> B 45 321
> B 47 512
> C 61 521
> C 68 235
> 
> Thank you
> 
> On Wed, Dec 6, 2017 at 10:34 PM, David Winsemius  
> wrote:
>> 
>>> On Dec 6, 2017, at 4:27 PM, Ashta  wrote:
>>> 
>>> Thank you Ista! Worked fine.
>> 
>> Here's another (possibly more direct in its logic?):
>> 
>> DM[ !ave(DM$x, DM$GR, FUN= function(x) {!length(unique(x))==1}), ]
>>  GR  x   y
>> 5  B 25 321
>> 6  B 25 512
>> 7  B 25 123
>> 8  B 25 451
>> 
>> --
>> David
>> 
>>> On Wed, Dec 6, 2017 at 5:59 PM, Ista Zahn  wrote:
>>>> Hi Ashta,
>>>> 
>>>> There are many ways to do it. Here is one:
>>>> 
>>>> vars <- sapply(split(DM$x, DM$GR), var)
>>>> DM[DM$GR %in% names(vars[vars > 0]), ]
>>>> 
>>>> Best
>>>> Ista
>>>> 
>>>> On Wed, Dec 6, 2017 at 6:58 PM, Ashta  wrote:
>>>>> Thank you Jeff,
>>>>> 
>>>>> subset( DM, "B" != x ), this works if I know the group only.
>>>>> But if I don't know that group in this case "B", how do I identify
>>>>> group(s) that  all elements of x have the same value?
>>>>> 
>>>>> On Wed, Dec 6, 2017 at 5:48 PM, Jeff Newmiller  
>>>>> wrote:
>>>>>> subset( DM, "B" != x )
>>>>>> 
>>>>>> This is covered in the Introduction to R document that comes with R.
>>>>>> --
>>>>>> Sent from my phone. Please excuse my brevity.
>>>>>> 
>>>>>> On December 6, 2017 3:21:12 PM PST, David Winsemius 
>>>>>>  wrote:
>>>>>>> 
>>>>>>>> On Dec 6, 2017, at 3:15 PM, Ashta  wrote:
>>>>>>>> 
>>>>>>>> Hi all,
>>>>>>>> In a data set I have group(GR) and two variables   x and y. I want to
>>>>>>>> remove a  group that have  the same record for the x variable in each
>>>>>>>> row.
>>>>>>>> 
>>>>>>>> DM <- read.table( text='GR x y
>>>>>>>> A 25 125
>>>>>>>> A 23 135
>>>>>>>> A 14 145
>>>>>>>> A 12 230
>>>>>>>> B 25 321
>>>>>>>> B 25 512
>>>>>>>> B 25 123
>>>>>>>> B 25 451
>>>>>>>> C 11 521
>>>>>>>> C 14 235
>>>>>>>> C 15 258
>>>>>>>> C 10 654',header = TRUE, stringsAsFactors = FALSE)
>>>>>>>> 
>>>>>>>> In this example the output should contain group A and C  as group B
>>>>>>>> has   the same record  for the variable x .
>>>>>>>> 
>>>>>>>> The result will be
>>>>>>>> A 25 125
>>>>>>>> A 23 135
>>>>>>>> A 14 145
>>>>>>>> A 12 230
>>>>>>>> C 11 521
>>>>>>>> C 14 235
>>>>>>>> C 15 258
>>>>>>>> C 10 654
>>>>>>> 
>>>>>>> Try:
>>>>>>> 
>>>>>>> DM[ !duplicated(DM$x) , ]
>>>>>>>> 
>>&g

Re: [R] Remove

[David Winsemius]

> On Dec 8, 2017, at 6:16 PM, David Winsemius  wrote:
> 
> 
>> On Dec 8, 2017, at 4:48 PM, Ashta  wrote:
>> 
>> Hi David, Ista and all,
>> 
>> I  have one related question  Within one group I want to keep records
>> conditionally.
>> example within
>> group A I want keep rows that have  " x" values  ranged  between 15 and 30.
>> group B I want keep rows that have  " x" values  ranged  between  40 and 50.
>> group C I want keep rows that have  " x" values  ranged  between  60 and 75.
> 
> When you have a problem where there are multiple "parallel: parameters, the 
> function to "reach for" is `mapply`. 
> 
>mapply( your_selection_func, group_vec, min_vec, max_vec)
> 
> ... and this will probably return the values as a list (of dataframes if you 
> build the function correctly,  so you may may need to then do:
> 
>do.call(rbind, ...)

 do.call( rbind, 
mapply( function(dat, grp, minx, maxx) {dat[ dat$GR==grp & dat$x >= minx & 
dat$x <= maxx, ]}, 
grp=LETTERS[1:3], minx=c(15,40,60), maxx=c(30,50,75) ,
MoreArgs=list(dat=DM),
IMPLIFY=FALSE))
 GR  x   y
A.1   A 25 125
A.2   A 23 135
B.5   B 45 321
B.6   B 47 512
C.9   C 61 521
C.10  C 68 235

> 
> -- 
> David.
>> 
>> 
>> DM <- read.table( text='GR x y
>> A 25 125
>> A 23 135
>> A 14 145
>> A 35 230
>> B 45 321
>> B 47 512
>> B 53 123
>> B 55 451
>> C 61 521
>> C 68 235
>> C 85 258
>> C 80 654',header = TRUE, stringsAsFactors = FALSE)
>> 
>> 
>> The end result will be
>> A 25 125
>> A 23 135
>> B 45 321
>> B 47 512
>> C 61 521
>> C 68 235
>> 
>> Thank you
>> 
>> On Wed, Dec 6, 2017 at 10:34 PM, David Winsemius  
>> wrote:
>>> 
>>>> On Dec 6, 2017, at 4:27 PM, Ashta  wrote:
>>>> 
>>>> Thank you Ista! Worked fine.
>>> 
>>> Here's another (possibly more direct in its logic?):
>>> 
>>> DM[ !ave(DM$x, DM$GR, FUN= function(x) {!length(unique(x))==1}), ]
>>> GR  x   y
>>> 5  B 25 321
>>> 6  B 25 512
>>> 7  B 25 123
>>> 8  B 25 451
>>> 
>>> --
>>> David
>>> 
>>>> On Wed, Dec 6, 2017 at 5:59 PM, Ista Zahn  wrote:
>>>>> Hi Ashta,
>>>>> 
>>>>> There are many ways to do it. Here is one:
>>>>> 
>>>>> vars <- sapply(split(DM$x, DM$GR), var)
>>>>> DM[DM$GR %in% names(vars[vars > 0]), ]
>>>>> 
>>>>> Best
>>>>> Ista
>>>>> 
>>>>> On Wed, Dec 6, 2017 at 6:58 PM, Ashta  wrote:
>>>>>> Thank you Jeff,
>>>>>> 
>>>>>> subset( DM, "B" != x ), this works if I know the group only.
>>>>>> But if I don't know that group in this case "B", how do I identify
>>>>>> group(s) that  all elements of x have the same value?
>>>>>> 
>>>>>> On Wed, Dec 6, 2017 at 5:48 PM, Jeff Newmiller 
>>>>>>  wrote:
>>>>>>> subset( DM, "B" != x )
>>>>>>> 
>>>>>>> This is covered in the Introduction to R document that comes with R.
>>>>>>> --
>>>>>>> Sent from my phone. Please excuse my brevity.
>>>>>>> 
>>>>>>> On December 6, 2017 3:21:12 PM PST, David Winsemius 
>>>>>>>  wrote:
>>>>>>>> 
>>>>>>>>> On Dec 6, 2017, at 3:15 PM, Ashta  wrote:
>>>>>>>>> 
>>>>>>>>> Hi all,
>>>>>>>>> In a data set I have group(GR) and two variables   x and y. I want to
>>>>>>>>> remove a  group that have  the same record for the x variable in each
>>>>>>>>> row.
>>>>>>>>> 
>>>>>>>>> DM <- read.table( text='GR x y
>>>>>>>>> A 25 125
>>>>>>>>> A 23 135
>>>>>>>>> A 14 145
>>>>>>>>> A 12 230
>>>>>>>>> B 25 321
>>>>>>>>> B 25 512
>>>>>>>>> B 25 123
>>>>>>>>> B 25 451
>>>>>>>>> C 11 521
>>>>>>>>> C 14 235
>>>>>>>>> C 15 258
&g

Re: [R] MSVAR model

[David Winsemius]

> On Dec 9, 2017, at 8:51 PM, ah a via R-help  wrote:
> 
> Hello,
> 
> As I'm interested to search about the monetary transmission channel in our
> country by MSVAR model,Could you do me favor and tell me How I can run
> different types of MSVAR model (such as MSIAH(2)-VAR(2)) and finding
> impulse response function in different regimes and also variance
> decomposition?

This is uncomfortably close to an existing fortune:

fortunes::fortune("brain surgery")

Why don't you read the Posting Guide. Then perhaps you can improve your 
question?

-- 
David.
> 
>   [[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.

David Winsemius
Alameda, CA, USA

'Any technology distinguishable from magic is insufficiently advanced.'   
-Gehm's Corollary to Clarke's Third Law

__
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] delta and sd parameters for power.t.test

[David Winsemius]

> On Dec 14, 2017, at 2:26 PM, Lasse Kliemann  wrote:
> 
> What is the rationale behind having both the delta and sd parameters for
> the power.t.test function?

One is the standard deviation of the hypothesized data (or pooled sd in the 
case of two sample) under the "alternative" and one is the mean of that data 
(or equivalently the differences if this is a two-sample test).

> For the relevant noncentrality parameter, we
> only need the ratio delta/sd. If my effect size is given as Cohen's d,
> then I only got that ratio and not sd.
> 

If the data has been "standardized", then Cohen's d could be given to the 
function as the value for delta since the sd default is 1.


> As far as I see, in such a case, I can specify delta=d and leave sd at
> its default value 1. Is this correct or am I missing something?
> 

Actually it sounds as though we are missing something. The power.t.test 
function makes no mention of Cohen's d or effect size. Is this question in 
response to a homework assignment about which we have not been informed?

Further questions should have some actual R code to make this an on-topic 
discussion for Rhelp.


> Thanks.
> __
> 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.

David Winsemius
Alameda, CA, USA

'Any technology distinguishable from magic is insufficiently advanced.'   
-Gehm's Corollary to Clarke's Third Law

__
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] something weird has happened....!!!!!!!!!!

[David Winsemius]

> On Dec 15, 2017, at 4:45 AM, akshay kulkarni  wrote:
> 
> dear Members,
> 
> 
> 
> Today something weird has happened on my R console. I have attached two 
> screenshots of the same vector in my R console but they differ.
> 
> 
> Also one of my function returns negative values, even after double checking 
> the code, which should return only positive values..
> 
> 
> Whats wrong..? Reinstall R?
> 
> 
> Thanks for help

We have no way of knowing what the function `yguii` might be doing. Apparently 
it has some sort of "random" internal process. If you want to see whether this 
is really a deterministic process (a pseudo-random process), then you should 
use `set.seed` to establish a reproducible state in the RNG prior to your 
function call.
> 
> 
> AKSHAY M KULKARNI
> 
> 
> __
> 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.

David Winsemius
Alameda, CA, USA

'Any technology distinguishable from magic is insufficiently advanced.'   
-Gehm's Corollary to Clarke's Third Law

__
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] Errors in reading in txt files

[David Winsemius]

> On Dec 15, 2017, at 9:21 AM, lily li  wrote:
> 
> I use the method, df$Time = as.POSIXct(df$Time), but it has the warning
> message:
> Error in as.POSIXlt.character(x, tz, ...) :
>  character string is not in a standard unambiguous format

That's because your date-time data is not in "%Y-%m-%d %H:%M" format. Read:

 ?strptime

-- 
David.
> 
> On Thu, Dec 14, 2017 at 1:31 PM, MacQueen, Don  wrote:
> 
>> In addition to which, I would recommend
>> 
>> df <- read.table("DATAM", header = TRUE, fill = TRUE,
>> stringsAsFactors=FALSE)
>> 
>> and then converting the Time column to POSIXct date-time values using
>>  as.POSIXct()
>> specifying the format using formatting codes found in
>>  ?strptime
>> because the times are not in the POSIXct default format.
>> 
>> 
>> This example might indicate the idea:
>> 
>>> as.POSIXct('2012-10-12 13:14')
>> [1] "2012-10-12 13:14:00 PDT"
>>> class(as.POSIXct('2012-10-12 13:14'))
>> [1] "POSIXct" "POSIXt"
>> 
>> -Don
>> 
>> --
>> Don MacQueen
>> Lawrence Livermore National Laboratory
>> 7000 East Ave., L-627
>> Livermore, CA 94550
>> 925-423-1062
>> Lab cell 925-724-7509
>> 
>> 
>> 
>> On 12/14/17, 11:01 AM, "R-help on behalf of Ista Zahn" <
>> r-help-boun...@r-project.org on behalf of istaz...@gmail.com> wrote:
>> 
>>On Thu, Dec 14, 2017 at 1:58 PM, Berend Hasselman 
>> wrote:
>>> 
>>>> On 14 Dec 2017, at 19:36, lily li  wrote:
>>>> 
>>>> Hi R users,
>>>> 
>>>> I have a question about reading from text files. The file has the
>> structure
>>>> below:
>>>> 
>>>> TimeColumn1   Column2
>>>> 01.01.2001-12:00:00
>>> 
>>> This line does not contain 3 elements; only one.
>>> You'll have to fix that line. Delete it, prepend it with a comment
>> character of add enough columns.
>> 
>>I definitely don't recommend that. Instead, read
>> 
>>?read.table
>> 
>>to learn about the "fill" and "header" arguments.
>> 
>>df = read.table("DATAM", header = TRUE, fill = TRUE)
>> 
>>will probably work.
>> 
>>Best,
>>Ista
>> 
>> 
>>> 
>>> 
>>> Berend
>>> 
>>>> 01.01.2001-24:00:0012 11
>>>> 01.02.2001-12:00:0013 10
>>>> 01.02.2001-24:00:0011 12
>>>> 01.03.2001-12:00:0015 11
>>>> 01.03.2001-24:00:0016 10
>>>> ...
>>>> 
>>>> I just use the simple script to open it: df = read.table('DATAM',
>> head=T).
>>>> 
>>>> But it has the error and thus cannot read the file:
>>>> Error in scan(file = file, what = what, sep = sep, quote = quote,
>> dec =
>>>> dec,  :
>>>> line 1 did not have 3 elements
>>>> 
>>>> How to read it with three fixed columns, and how to read the time
>> format in
>>>> the first column correctly? Thanks for your help.
>>>> 
>>>>  [[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.
>> 
>> 
>> 
> 
>   [[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.

David Winsemius
Alameda, CA, USA

'Any technology distinguishable from magic is insufficiently advanced.'   
-Gehm's Corollary to Clarke's Third Law

__
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] offset with a factor

[David Winsemius]

> On Dec 20, 2017, at 6:51 AM, Bond, Stephen  wrote:
> 
> Knowledgeable useRs,
> 
> Please, advise how to use offset with a factor. I estimate monthly effects 
> from a much bigger data set as monthly effects seem to be stable, and other 
> variables are estimated from a small, but recent data set as there is 
> variation in those non-seasonal coefficients.
> How can I use the seasonality estimates from the big data set as an offset 
> provided to the small data set. I know an offset is supposed to be 
> quantitative, but this is such a practical and sensible scenario, I feel 
> compelled. Assume I have 11 coefs estimated with contr.sum.

Why not create a variable that specifies the relevant coefficient for each of 
the various levels of the factor (and 0 for the reference level? Then pass that 
variable to the offset argument.

-- 
David.
> 
> Thanks everybody
> 
> Stephen
> 
> 
>   [[alternative HTML version deleted]]

Do note: Rhelp is a plain-text mailing list.
> 
> __
> 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.

David Winsemius
Alameda, CA, USA

'Any technology distinguishable from magic is insufficiently advanced.'   
-Gehm's Corollary to Clarke's Third Law

__
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] problem in installing "simpleaffy"

[David Winsemius]

> On Dec 20, 2017, at 7:55 AM, Rahele Amirkhah  
> wrote:
> 
> Dear Madam/ Sir,
> I am using R version 3.4.2. I want to analyse microarray data. when I want to 
> install "simpleaffy" package I get this error "package ‘simpleaffy’ is not 
> available (for R version 3.4.2)". I have the same problem with R version 
> 3.3.2.

That's because it's not a CRAN package. It's a Bioc package.


> Could you please help me to solve it?
> I am working with RStudio 0.99.903.exe. I also have problem in getting the 
> new release of RStudio. Could you please tell me how I can get the last 
> release? 

That's off-topic for this mailing list. Ask at RStudio

> Many thanks in advance.Best regards,Raheleh Amirkhah
>   [[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.

David Winsemius
Alameda, CA, USA

'Any technology distinguishable from magic is insufficiently advanced.'   
-Gehm's Corollary to Clarke's Third Law

__
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] LSD-test

[David Winsemius]

> On Dec 25, 2017, at 10:21 AM, Ahmed Attia  wrote:
> 
> LSD.test

?LSD.test
No documentation for ‘LSD.test’ in specified packages and libraries:
you could try ‘??LSD.test’

-- 
David Winsemius
Alameda, CA, USA

'Any technology distinguishable from magic is insufficiently advanced.'   
-Gehm's Corollary to Clarke's Third Law

__
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] LSD-test

[David Winsemius]

> On Dec 25, 2017, at 2:09 PM, Ahmed Attia  wrote:
> 
> The model should be class aov or lm and my model class is aovlist.
> tried tidy from broom library but did not work. To make it class aov,
> I had to remove the error term;
> 
> model <- 
> aov(Rotationdata_R$`GY(Mg/ha)`~Rep+code*as.factor(Nitrogen),data=Rotationdata_R)

You seemed to have missed my point that LSD.test is not in the packages loaded 
by default.


> Ahmed Attia, Ph.D.
> Agronomist & Soil Scientist
> 
> 
> 
> 
> 
> 
> On Mon, Dec 25, 2017 at 7:38 PM, David Winsemius  
> wrote:
>> 
>>> On Dec 25, 2017, at 10:21 AM, Ahmed Attia  wrote:
>>> 
>>> LSD.test
>> 
>> ?LSD.test
>> No documentation for ‘LSD.test’ in specified packages and libraries:
>> you could try ‘??LSD.test’
>> 
>> --
>> David Winsemius
>> Alameda, CA, USA
>> 
>> 'Any technology distinguishable from magic is insufficiently advanced.'   
>> -Gehm's Corollary to Clarke's Third Law
>> 
>> 
>> 
>> 
>> 

David Winsemius
Alameda, CA, USA

'Any technology distinguishable from magic is insufficiently advanced.'   
-Gehm's Corollary to Clarke's Third Law

__
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] RQuantLib

[David Winsemius]

> On Dec 30, 2017, at 7:54 AM, rsherry8  wrote:
> 
> OA,
> 
> Thanks for the response. I downloaded the file RQuantLib_0.4.4.tar.gz 
> into a directory c:\r.zip on my Windows machine. I then ran the 
> following command and I received the following output:
> 
>> install.packages("RQuantLib", lib="/r.zip/")
> Warning message:
> package ‘RQuantLib’ is not available (for R version 3.4.3)

The install.packages command should include repo=NULL when installing from 
local binary package, and it also should include type="source when the package 
is not binary.
> 
> I did not unpack the .gz file. Should I have?
> 
> Please comment.
> Bob Sherry
> 
> On 12/30/2017 2:24 AM, Orvalho Augusto wrote:
>> Hi Bob,
>> 
>> I don't know what is the cause of your trouble but try this:
>> 1. Download the zip of package.
>> 
>> 2. And install it from local zip files. This you find on the Packages 
>> menu.
>> 
>> Hope it helps
>> OA
>> 
>> 
>> On Fri, Dec 29, 2017 at 10:31 AM, rsherry8 > <mailto:rsher...@comcast.net>> wrote:
>> 
>>Joshua,
>> 
>>Thanks for the response. When you said at least version 3.4.0, I
>>upgraded to 3.4.2 which I believe is the current version. Now, I
>>attempted to install the package RQuantLib but it did not work.
>>Here is what I got:
>> 
>>> install.packages("RQuantLib")
>>Installing package into ‘C:/Users/rsher/Documents/R/win-library/3.4’
>>(as ‘lib’ is unspecified)
>>Warning message:
>>package ‘RQuantLib’ is not available (for R version 3.4.2)
>> 
>> 
>>Please help.
>>Thanks,
>>Bob Sherry
>> 
>> 
>>On 12/28/2017 10:28 PM, Joshua Ulrich wrote:
>> 
>>On Thu, Dec 28, 2017 at 6:02 PM, rsherry8
>>mailto:rsher...@comcast.net>> wrote:
>> 
>>I have recently installed R on my new computer. I also
>>want to install the
>>package RQuantLib. So I run the following command and get
>>the following
>>output:
>> 
>>  install.packages("RQuantLib")
>> 
>>Installing package into
>>‘C:/Users/rsher/Documents/R/win-library/3.2’
>>(as ‘lib’ is unspecified)
>>--- Please select a CRAN mirror for use in this session ---
>>Warning message:
>>package ‘RQuantLib’ is not available (for R version 3.2.4
>>Revised)
>> 
>>The package did not install. Am I doing something wrong.
>>Is the package
>>going to be updated for the latest version of R?
>> 
>>Windows binary packages are only built for the most current
>>(major)
>>version of R.  You need to upgrade to at least R-3.4.0, or you
>>will
>>have to install RQuantLib (and therefore QuantLib itself) from
>>source.
>> 
>>Thanks,
>>Bob
>> 
>>__
>>R-help@r-project.org <mailto:R-help@r-project.org> mailing
>>list -- To UNSUBSCRIBE and more, see
>>https://stat.ethz.ch/mailman/listinfo/r-help
>><https://stat.ethz.ch/mailman/listinfo/r-help>
>>PLEASE do read the posting guide
>>http://www.R-project.org/posting-guide.html
>><http://www.R-project.org/posting-guide.html>
>>and provide commented, minimal, self-contained,
>>reproducible code.
>> 
>> 
>> 
>> 
>>__
>>R-help@r-project.org <mailto:R-help@r-project.org> mailing list --
>>To UNSUBSCRIBE and more, see
>>https://stat.ethz.ch/mailman/listinfo/r-help
>><https://stat.ethz.ch/mailman/listinfo/r-help>
>>PLEASE do read the posting guide
>>http://www.R-project.org/posting-guide.html
>><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.

David Winsemius
Alameda, CA, USA

'Any technology distinguishable from magic is insufficiently advanced.'   
-Gehm's Corollary to Clarke's Third Law

__
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] RQuantLib

[David Winsemius]

> On Dec 30, 2017, at 9:44 AM, David Winsemius  wrote:
> 
> 
>> On Dec 30, 2017, at 7:54 AM, rsherry8  wrote:
>> 
>> OA,
>> 
>> Thanks for the response. I downloaded the file RQuantLib_0.4.4.tar.gz 
>> into a directory c:\r.zip on my Windows machine. I then ran the 
>> following command and I received the following output:
>> 
>>> install.packages("RQuantLib", lib="/r.zip/")
>> Warning message:
>> package ‘RQuantLib’ is not available (for R version 3.4.3)
> 
> The install.packages command should include repo=NULL when installing from 
> local binary package, and it also should include type="source when the 
> package is not binary.

Furthermore, the system requirements in the DESCRIPTION file are:


SystemRequirements: QuantLib library (>= 1.8.0) from http://quantlib.org, 
Boost library from http://www.boost.org


>> 
>> I did not unpack the .gz file. Should I have?
>> 
>> Please comment.
>> Bob Sherry
>> 
>> On 12/30/2017 2:24 AM, Orvalho Augusto wrote:
>>> Hi Bob,
>>> 
>>> I don't know what is the cause of your trouble but try this:
>>> 1. Download the zip of package.
>>> 
>>> 2. And install it from local zip files. This you find on the Packages 
>>> menu.
>>> 
>>> Hope it helps
>>> OA
>>> 
>>> 
>>> On Fri, Dec 29, 2017 at 10:31 AM, rsherry8 >> <mailto:rsher...@comcast.net>> wrote:
>>> 
>>>   Joshua,
>>> 
>>>   Thanks for the response. When you said at least version 3.4.0, I
>>>   upgraded to 3.4.2 which I believe is the current version. Now, I
>>>   attempted to install the package RQuantLib but it did not work.
>>>   Here is what I got:
>>> 
>>>> install.packages("RQuantLib")
>>>   Installing package into ‘C:/Users/rsher/Documents/R/win-library/3.4’
>>>   (as ‘lib’ is unspecified)
>>>   Warning message:
>>>   package ‘RQuantLib’ is not available (for R version 3.4.2)
>>> 
>>> 
>>>   Please help.
>>>   Thanks,
>>>   Bob Sherry
>>> 
>>> 
>>>   On 12/28/2017 10:28 PM, Joshua Ulrich wrote:
>>> 
>>>   On Thu, Dec 28, 2017 at 6:02 PM, rsherry8
>>>   mailto:rsher...@comcast.net>> wrote:
>>> 
>>>   I have recently installed R on my new computer. I also
>>>   want to install the
>>>   package RQuantLib. So I run the following command and get
>>>   the following
>>>   output:
>>> 
>>> install.packages("RQuantLib")
>>> 
>>>   Installing package into
>>>   ‘C:/Users/rsher/Documents/R/win-library/3.2’
>>>   (as ‘lib’ is unspecified)
>>>   --- Please select a CRAN mirror for use in this session ---
>>>   Warning message:
>>>   package ‘RQuantLib’ is not available (for R version 3.2.4
>>>   Revised)
>>> 
>>>   The package did not install. Am I doing something wrong.
>>>   Is the package
>>>   going to be updated for the latest version of R?
>>> 
>>>   Windows binary packages are only built for the most current
>>>   (major)
>>>   version of R.  You need to upgrade to at least R-3.4.0, or you
>>>   will
>>>   have to install RQuantLib (and therefore QuantLib itself) from
>>>   source.
>>> 
>>>   Thanks,
>>>   Bob
>>> 
>>>   __
>>>   R-help@r-project.org <mailto:R-help@r-project.org> mailing
>>>   list -- To UNSUBSCRIBE and more, see
>>>   https://stat.ethz.ch/mailman/listinfo/r-help
>>>   <https://stat.ethz.ch/mailman/listinfo/r-help>
>>>   PLEASE do read the posting guide
>>>   http://www.R-project.org/posting-guide.html
>>>   <http://www.R-project.org/posting-guide.html>
>>>   and provide commented, minimal, self-contained,
>>>   reproducible code.
>>> 
>>> 
>>> 
>>> 
>>>   __
>>>   R-help@r-project.org <mailto:R-help@r-project.org> mailing list --
>>>   To UNSUBSCRIBE and more, see
>>>   https://stat.ethz.ch/mailman/listinfo/r-help
>>>   <https://stat.ethz.ch/

Re: [R] RQuantLib

[David Winsemius]

> On Dec 30, 2017, at 9:46 AM, David Winsemius  wrote:
> 
> 
>> On Dec 30, 2017, at 9:44 AM, David Winsemius  wrote:
>> 
>> 
>>> On Dec 30, 2017, at 7:54 AM, rsherry8  wrote:
>>> 
>>> OA,
>>> 
>>> Thanks for the response. I downloaded the file RQuantLib_0.4.4.tar.gz 
>>> into a directory c:\r.zip on my Windows machine. I then ran the 
>>> following command and I received the following output:
>>> 
>>>> install.packages("RQuantLib", lib="/r.zip/")
>>> Warning message:
>>> package ‘RQuantLib’ is not available (for R version 3.4.3)
>> 
>> The install.packages command should include repo=NULL when installing from 
>> local binary package, and it also should include type="source"

Added the missing dbl-quote.

>>  when the package is not binary.
> 
> Furthermore, the system requirements in the DESCRIPTION file are:
> 
> 
> SystemRequirements:   QuantLib library (>= 1.8.0) from http://quantlib.org, 
> Boost library from http://www.boost.org

And finally (perhaps):

The binaries at CRAN for windows are still at 0.4.2 so you might consider 
install that version from source:

RQuantLib_0.4.2.tar.gz at 
https://cran.r-project.org/src/contrib/Archive/RQuantLib/RQuantLib_0.4.2.tar.gz


> 
> 
>>> 
>>> I did not unpack the .gz file. Should I have?
>>> 
>>> Please comment.
>>> Bob Sherry
>>> 
>>> On 12/30/2017 2:24 AM, Orvalho Augusto wrote:
>>>> Hi Bob,
>>>> 
>>>> I don't know what is the cause of your trouble but try this:
>>>> 1. Download the zip of package.
>>>> 
>>>> 2. And install it from local zip files. This you find on the Packages 
>>>> menu.
>>>> 
>>>> Hope it helps
>>>> OA
>>>> 
>>>> 
>>>> On Fri, Dec 29, 2017 at 10:31 AM, rsherry8 >>> <mailto:rsher...@comcast.net>> wrote:
>>>> 
>>>>  Joshua,
>>>> 
>>>>  Thanks for the response. When you said at least version 3.4.0, I
>>>>  upgraded to 3.4.2 which I believe is the current version. Now, I
>>>>  attempted to install the package RQuantLib but it did not work.
>>>>  Here is what I got:
>>>> 
>>>>> install.packages("RQuantLib")
>>>>  Installing package into ‘C:/Users/rsher/Documents/R/win-library/3.4’
>>>>  (as ‘lib’ is unspecified)
>>>>  Warning message:
>>>>  package ‘RQuantLib’ is not available (for R version 3.4.2)
>>>> 
>>>> 
>>>>  Please help.
>>>>  Thanks,
>>>>  Bob Sherry
>>>> 
>>>> 
>>>>  On 12/28/2017 10:28 PM, Joshua Ulrich wrote:
>>>> 
>>>>  On Thu, Dec 28, 2017 at 6:02 PM, rsherry8
>>>>  mailto:rsher...@comcast.net>> wrote:
>>>> 
>>>>  I have recently installed R on my new computer. I also
>>>>  want to install the
>>>>  package RQuantLib. So I run the following command and get
>>>>  the following
>>>>  output:
>>>> 
>>>>install.packages("RQuantLib")
>>>> 
>>>>  Installing package into
>>>>  ‘C:/Users/rsher/Documents/R/win-library/3.2’
>>>>  (as ‘lib’ is unspecified)
>>>>  --- Please select a CRAN mirror for use in this session ---
>>>>  Warning message:
>>>>  package ‘RQuantLib’ is not available (for R version 3.2.4
>>>>  Revised)
>>>> 
>>>>  The package did not install. Am I doing something wrong.
>>>>  Is the package
>>>>  going to be updated for the latest version of R?
>>>> 
>>>>  Windows binary packages are only built for the most current
>>>>  (major)
>>>>  version of R.  You need to upgrade to at least R-3.4.0, or you
>>>>  will
>>>>  have to install RQuantLib (and therefore QuantLib itself) from
>>>>  source.
>>>> 
>>>>  Thanks,
>>>>  Bob
>>>> 
>>>>  __
>>>>  R-help@r-project.org <mailto:R-help@r-project.org> mailing
>>>>  list -- To UNSUBSCRIBE and more, see
>>>>  https://stat.ethz.ch/mailman/listinfo/r-help
&g

Re: [R] httr::content without message

[David Winsemius]

> On Jan 2, 2018, at 9:30 AM, Roy Mendelssohn - NOAA Federal 
>  wrote:
> 
> Hi All:
> 
> I am using httr to download files form a service, in this case a .csv file.  
> When I use httr::content on the result,  I get a message.  Since this will be 
> in a package.  I want to suppress the message,  but haven't figured out how 
> to do so.
> 
> The following should reproduce the result:
> 
> myURL <- 
> 'https://coastwatch.pfeg.noaa.gov/erddap/griddap/erdMH1sstdmday.csvp?time[0:1:last]'
> r1 <- httr::GET(myURL)
> junk <- httr::content(r1)

Instead try:

junk <- suppressMessages(httr::content(r1))

> 

> when the last command is run, you get:
> 
> Parsed with column specification:
> cols(
>  `time (UTC)` = col_datetime(format = "")
> )
> 
> I want to suppress that output.
> 
> Thanks,
> 
> -Roy
> 
> **
> "The contents of this message do not reflect any position of the U.S. 
> Government or NOAA."
> **
> Roy Mendelssohn
> Supervisory Operations Research Analyst
> NOAA/NMFS
> Environmental Research Division
> Southwest Fisheries Science Center
> ***Note new street address***
> 110 McAllister Way
> Santa Cruz, CA 95060
> Phone: (831)-420-3666
> Fax: (831) 420-3980
> e-mail: roy.mendelss...@noaa.gov www: http://www.pfeg.noaa.gov/
> 
> "Old age and treachery will overcome youth and skill."
> "From those who have been given much, much will be expected" 
> "the arc of the moral universe is long, but it bends toward justice" -MLK Jr.
> 
> __
> 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.

David Winsemius
Alameda, CA, USA

'Any technology distinguishable from magic is insufficiently advanced.'   
-Gehm's Corollary to Clarke's Third Law

__
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] summary.rms help

[David Winsemius]

> On Jan 3, 2018, at 11:57 AM, Andras Farkas via R-help  
> wrote:
> 
> Dear All,
> using the example from the help of summary.rms
> 
> library(rms)
> n <- 1000# define sample size 
> set.seed(17) # so can reproduce the results 
> age<- rnorm(n, 50, 10) 
> blood.pressure <- rnorm(n, 120, 15) 
> cholesterol<- rnorm(n, 200, 25) 
> sex<- factor(sample(c('female','male'), n,TRUE)) 
> label(age)<- 'Age'  # label is in Hmisc 
> label(cholesterol)<- 'Total Cholesterol' 
> label(blood.pressure) <- 'Systolic Blood Pressure' 
> label(sex)<- 'Sex' 
> units(cholesterol)<- 'mg/dl'   # uses units.default in Hmisc 
> units(blood.pressure) <- 'mmHg' 
> # Specify population model for log odds that Y=1 
> L <- .4*(sex=='male') + .045*(age-50) + 
> (log(cholesterol - 10)-5.2)*(-2*(sex=='female') + 2*(sex=='male')) 
> # Simulate binary y to have Prob(y=1) = 1/[1+exp(-L)] 
> y <- ifelse(runif(n) < plogis(L), 1, 0) 
> ddist <- datadist(age, blood.pressure, cholesterol, sex) 
> options(datadist='ddist') 
> fit <- lrm(y ~ blood.pressure + sex * (age + rcs(cholesterol,4)))
> s <- summary(fit) 
> plot(s)
> as you will see the plot will by default include the low and high values from 
> the summary printed on the plot to the right of the variable name... Any 
> thoughts on how printing these low and high values can be suppressed, ie: 
> prevent them from being printed?
> 

Luke, ... Look at the code!

The values are suppressed if the "Diff." has NA's so ...

s[ , "Diff."] <- NA
plot(s)

-- 
David Winsemius
Alameda, CA, USA

'Any technology distinguishable from magic is insufficiently advanced.'   
-Gehm's Corollary to Clarke's Third Law

__
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] svm

[David Winsemius]

> On Jan 10, 2018, at 10:50 AM, AbouEl-Makarim Aboueissa 
>  wrote:
> 
> Dear All:
> 
> 
> I am trying to use the R function "svm" with  "type =C-classification" ,
> but I got the following error message
> 
> 
> SVM.Model1 <- svm(type ~ ., data=my.data.x1x2y, *type='C-classification'*,
> kernel='linear',scale=FALSE)
> 
> *Error in eval(predvars, data, env) : object 'type' not found*

Yopu misspelled the name of your "typey" variable.
> 
> 
> I am wondering if I should install a specific R package(s).
> 
> 
> 
> *Here is my codes:*
> 
> 
> feature.x1 <- c(0.25,0.5,1,1,1.5,2,2.25,2.5,2,1,3,   5,3.75,
> 1,3.5,4,4,5,5.5,6,6,6.5)
> 
> length(feature.x1)
> 
> 
> 
> feature.x2 <- c(2,3.5,1,2.5,1.75,2,1.75,1.5,2.5,1,1,
> 3.5,3.5,5.8,3,4,4.5,5,4,1,4,3)
> 
> length(feature.x2)
> 
> 
> y <- c(rep(-1,11), rep(1,11))
> 
> typey<-as.factor(y)
> 
> 
> my.data.x1x2y <- data.frame(feature.x1, feature.x2, typey)
> 
> my.data.x1x2y
> 
> 
> 
> install.packages("e1071")
> 
> library(e1071)
> 
> 
> 
> SVM.Model1 <- svm(type ~ ., data=my.data.x1x2y, type='C-classification',
> kernel='linear',scale=FALSE)
> 
> plot(my.data.x1x2y[,-3],col=(typey+3)/2, pch=18, xlim=c(-1,6), ylim=c(-1,6))
> 
> box(lwd = 2, col="darkgreen")
> 
> 
> 
> 
> 
> 
> 
> with many thanks
> abou
> __
> 
> 
> *AbouEl-Makarim Aboueissa, PhD*
> 
> *Professor of Statistics*
> 
> *Department of Mathematics and Statistics*
> *University of Southern Maine*
> 
>   [[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.

David Winsemius
Alameda, CA, USA

'Any technology distinguishable from magic is insufficiently advanced.'   
-Gehm's Corollary to Clarke's Third Law

__
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] Help with packages (methods, stats, stats4)

[David Winsemius]

gt;>>> 
>>>>>> (i want to install those packages in order to use the "DVstats" package)
>>>>>> 
>>>>>> (i have the latest version of R (3.4.3 ) and Rstudio(1.2.240) )
>>>>>> 
>>>>>> thank you
>>>>>> Regards
>>>>>> 
>>>>>> 
>>>>>  [[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/posti
>>>> ng-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.

David Winsemius
Alameda, CA, USA

'Any technology distinguishable from magic is insufficiently advanced.'   
-Gehm's Corollary to Clarke's Third Law

__
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] Time-dependent coefficients in a Cox model with categorical variants

[David Winsemius]

> On Jan 15, 2018, at 12:58 AM, Max Shell  wrote:
> 
> Suppose I have a dataset contain three variants, looks like
>> head(dta)
> 
>  SextumorsizeHistology   time status
>01.52  12.1000 0
>11.81  38.4000 0
> .
> 
> Sex: 1 for male; 0 for female., two levels
> Histology: 1 for SqCC; 2 for High risk AC; 3 for low risk AC, three levels
> Now I need to get a Time-dependent coefficients cox fit:
> 
> library(survival)
> for(i in c(1,3) dta[,i] <- factor(dta[,i])
> fit <-
>  coxph(
>Surv(time, status) ~  Sex + tumorsize +  Histology + tt(Histology),
>data = dta,
>tt = function(x, t, ...) x * log(t)
>  )
> 
> But  I keep gettting this error says:
> 
> Error in if (any(infs)) warning(paste("Loglik converged before variable ",  :
>  missing value where TRUE/FALSE needed
> In addition: Warning message:
> In Ops.factor(x, log(t)) : ‘*’ not meaningful for factors.

The error message seems pretty clear. You are passing a factor to the x 
parameter of the tt function, and then you are attempting to multiply that 
value times log(t). You are lucky that it was a factor and not jsut an integer 
because then you might not have gotten an error or a warning.  I worry that 
your hopes of separate estimates may not be easily supportable by the 
tt-mechanism. However, the example of its use in the help page for `coxph` 
shows a spline function being passed and the boundary knots and weights mush 
estimated, so my fears may be over-blown. The knot locations and weights are 
not reported in the print.coxph but it doesn't look too difficult to extract 
then from the attributes of the model. 

You might look at the mechanism for estimation of spline component effects to 
see if you can learn how to estimate multiple components:

> coef( coxph(Surv(time, status) ~ ph.ecog + tt(age), data=lung,
+  tt=function(x,t,...) pspline(x + t/365.25))  )
   ph.ecog  ps(x + t/365.25)3  ps(x + t/365.25)4  ps(x + t/365.25)5  
ps(x + t/365.25)6  ps(x + t/365.25)7 
 0.4528363  0.2426635  0.4876185  0.7796924 
 1.0160954  1.0765967 
 ps(x + t/365.25)8  ps(x + t/365.25)9 ps(x + t/365.25)10 ps(x + t/365.25)11 
ps(x + t/365.25)12 ps(x + t/365.25)13 
 1.0449439  0.9170725  0.9276695  1.1349794 
 1.2837341  1.7024045 
ps(x + t/365.25)14 
 2.1863712 

> attr( coxph(Surv(time, status) ~ ph.ecog + tt(age), data=lung,
+  tt=function(x,t,...) pspline(x + t/365.25))$terms, "predvars" )
list(Surv(time, status), ph.ecog, pspline(age, nterm = 10, intercept = FALSE, 
Boundary.knots = c(39.0136892539357, 82.0848733744011), `NA` = NULL))

The coefficients for the spline also get reported as exponentiated values in 
the summary output. And if you  used a crossing operator in the formula you get 
some sort of interaction result. Whether it has any sensible interpretation is 
decision that's above my pay grade.

The code for `pspline` is readily available. It's not even necessary to sue the 
triple-colon or getAnywhere functions.

-- 
David.
> 
> How can I fix it? I know that the "Sex" and "Histology" are both
> categorical variants. I  want to have a model that have two β(t) = a +
> blog(t) for each histology level.
> Thank you！
> 
> __
> 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.

David Winsemius
Alameda, CA, USA

'Any technology distinguishable from magic is insufficiently advanced.'   
-Gehm's Corollary to Clarke's Third Law

__
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] reading lisp file in R

[David Winsemius]

> On Jan 17, 2018, at 8:22 PM, Ranjan Maitra  wrote:
> 
> Dear friends,
> 
> Is there a way to read data files written in lisp into R? 
> 
> Here is the file: 
> https://archive.ics.uci.edu/ml/machine-learning-databases/university/university.data
> 
> I would like to read it into R. Any suggestions?

It's just a text file. What difficulties are you having?
>  
> 
> Thanks very much in advance for pointers on this and best wishes,
> Ranjan
> 
> -- 
> Important Notice: This mailbox is ignored: e-mails are set to be deleted on 
> receipt. Please respond to the mailing list if appropriate. For those needing 
> to send personal or professional e-mail, please use appropriate addresses.
> 
> __
> 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.

David Winsemius
Alameda, CA, USA

'Any technology distinguishable from magic is insufficiently advanced.'   
-Gehm's Corollary to Clarke's Third Law

__
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 search r-help?

[David Winsemius]

> On Jan 19, 2018, at 10:57 AM, Mark Dwyer via R-help  
> wrote:
> 
> Also https://www.r-project.org/posting-guide.html indicates that 
> RSiteSearch() within R searches R-help but in my install (3.4.3) 
> RSiteSearch() only searches  "help pages, vignettes or taskviews"
> 
> 
> On 19/01/18 10:47, Mark Dwyer wrote:
>> 
>> I am new to this listand am unable to get the search tools listed on 
>> https://stat.ethz.ch/mailman/listinfo/r-help towork. What do people use to 
>> search the help archives?
>> 
>> 1. The google search box on http://tolstoy.newcastle.edu.au/~rking/R/
>>returns a 404 error.
>> 2. The http://finzi.psych.upenn.edu/ site has many references but I
>>don't see how to search r-help from there.
>> 3. The http://www.mail-archive.com/r-help@stat.math.ethz.ch/ site
>>seems to stop back in December 2007.
>> 4. I cannot get anything useful from
>>http://dir.gmane.org/gmane.comp.lang.r.general.

For R-help I've been using:

http://markmail.org/search/?q=list%3Aorg.r-project.r-help

To find documentation, I use:

install.packages("sos")
library(sos)
findFn("topic items")

You could conceivably use Googles advanced search with the site set to the 
"true" Archives.

https://www.google.com/search?as_q=words+to+search&as_epq=&as_oq=&as_eq=&as_nlo=&as_nhi=&lr=&cr=&as_qdr=all&as_sitesearch=https%3A%2F%2Fstat.ethz.ch%2Fpipermail%2Fr-help%2F&as_occt=any&safe=images&as_filetype=&as_rights=

-- 
David Winsemius
Alameda, CA, USA

'Any technology distinguishable from magic is insufficiently advanced.'   
-Gehm's Corollary to Clarke's Third Law

__
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] Help in Plotting in "fArma" Package

[David Winsemius]

The documentation say that additional arguments will be passed. I suspect this 
will be a base graphics plot. You should look at the code of plot.rsfit to 
determine which arguments get processed. 

Sent from my iPhone

> On Jan 25, 2018, at 10:30 AM, Moyukh Laha  wrote:
> 
> Hello,
>I am new to R and for some of my research work I am using 'fArma'
> package to estimate the Hurst parameter of a time series.
>When I am ding the following command :
> rsFit(data, doplot = TRUE)
> I am getting the R/S plot for that time series with default plot title,
> font size. However, I want to  change the axis size, font size etc of this
> plot, which I am unable to do as there is no formal argument here. Can
> anyone suggest something about this? How can I change the font size, axis
> size etc here?
> Thanks.
> 
>[[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] Help in Plotting in "fArma" Package

[David Winsemius]

> On Jan 26, 2018, at 9:51 AM, MacQueen, Don  wrote:
> 
> What Dave said, plus here's a hint. Try this example (which uses base 
> graphics):
> 
> plot(1:5)
> plot(1:5, cex.lab=2)
> 
> Then look at the help page for par
>   help('par')
> or
>   ?par
> to search for other graphics parameters (base graphics) you can use to change 
> various things.
> 
> Success will depend, as Dave indicated, on how the package author handled the 
> plotting options in rsFit().

Actually it's an S4 method:

showMethods("show", classes="fHURST", includeDefs=TRUE)
Function: show (package methods)
object="fHURST"
function (object) 
{
x = object
doplot = TRUE
cat("\nTitle:\n ", x@title, "\n", sep = "")
cat("\nCall:\n ")
cat(paste(deparse(x@call), sep = "\n", collapse = "\n"), 
"\n", sep = "")
... snipped

IAnd not easily susceptible to throwing base graphics parameters at ti since 
it's all hard-coded. I've suggested to the OP that hacking the show method is 
an accessible option.

-- 
David.

> 
> -Don
> 
> --
> Don MacQueen
> Lawrence Livermore National Laboratory
> 7000 East Ave., L-627
> Livermore, CA 94550
> 925-423-1062
> Lab cell 925-724-7509
> 
> 
> On 1/25/18, 5:56 PM, "R-help on behalf of David Winsemius" 
>  wrote:
> 
>The documentation say that additional arguments will be passed. I suspect 
> this will be a base graphics plot. You should look at the code of plot.rsfit 
> to determine which arguments get processed. 
> 
>Sent from my iPhone
> 
>> On Jan 25, 2018, at 10:30 AM, Moyukh Laha  wrote:
>> 
>> Hello,
>>   I am new to R and for some of my research work I am using 'fArma'
>> package to estimate the Hurst parameter of a time series.
>>   When I am ding the following command :
>>rsFit(data, doplot = TRUE)
>> I am getting the R/S plot for that time series with default plot title,
>> font size. However, I want to  change the axis size, font size etc of this
>> plot, which I am unable to do as there is no formal argument here. Can
>> anyone suggest something about this? How can I change the font size, axis
>> size etc here?
>> Thanks.
>> 
>>   [[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.
> 
> 

David Winsemius
Alameda, CA, USA

'Any technology distinguishable from magic is insufficiently advanced.'   
-Gehm's Corollary to Clarke's Third Law

__
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] Fortune candidate

[David Winsemius]

John (to a serial querulant):  

...but with such a sweeping lack of
information from you, don't congratulate yourself if you get a helpful
answer.  It wasn't your fault.


David Winsemius
Alameda, CA, USA

'Any technology distinguishable from magic is insufficiently advanced.'   
-Gehm's Corollary to Clarke's Third Law

__
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] Newbie wants to compare 2 huge RDSs row by row.

[David Winsemius]

> On Jan 27, 2018, at 1:18 PM, Marsh Hardy ARA/RISK  wrote:
> 
> Hi Guys, I apologize for my rank & utter newness at R.
> 
> I used summary() and found about 95 variables, both character and numeric, 
> all with "Length:368842" I assume is the # of records.
> 
> I'd like to know the record number (row #?) of any record where the data 
> doesn't match in the 2 files of what should be the same output.

The 'length' function returns the number of items at the top level of a list. 
As such it returns the number of columns of a dataframe which is a type of 
list. The quotes around "Length:368842" make it appear that you are quoting 
from some sort of output, although its nature is unclear. Rather than 
abstracting from information that you don't know how to interpret it would have 
been better to include the entire output. The `summary` function is generic and 
therefore can vary widely in what it prints to the console.

You should instead post the output from str() applied to each of your "files". 
It does a better job of displaying object structure.

-- 
David.


> Thanks in advance, M.
> 
> //
> 
> From: Ulrik Stervbo [ulrik.ster...@gmail.com]
> Sent: Saturday, January 27, 2018 10:00 AM
> To: Eric Berger
> Cc: Marsh Hardy ARA/RISK; r-help@r-project.org
> Subject: Re: [R] Newbie wants to compare 2 huge RDSs row by row.
> 
> Also, it will be easier to provide helpful information if you'd describe what 
> in your data you want to compare and what you hope to get out of the 
> comparison.
> 
> Best wishes,
> Ulrik
> 
> Eric Berger mailto:ericjber...@gmail.com>> schrieb am 
> Sa., 27. Jan. 2018, 08:18:
> Hi Marsh,
> An RDS is not a data structure such as a data.frame. It can be anything.
> For example if I want to save my objects a, b, c I could do:
>> saveRDS( list(a,b,c,), file="tmp.RDS")
> Then read them back later with
>> myList <- readRDS( "tmp.RDS" )
> 
> Do you have additional information about your "RDSs" ?
> 
> Eric
> 
> 
> On Sat, Jan 27, 2018 at 6:54 AM, Marsh Hardy ARA/RISK 
> mailto:mha...@ara.com>>
> wrote:
> 
>> Each RDS is 40 MBs. What's a slick code to compare them row by row, IDing
>> row numbers with mismatches?
>> 
>> Thanks in advance.
>> 
>> //
>> 
>> __
>> R-help@r-project.org<mailto: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<mailto: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.

David Winsemius
Alameda, CA, USA

'Any technology distinguishable from magic is insufficiently advanced.'   
-Gehm's Corollary to Clarke's Third Law

__
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] semPLS package will not load seems to be failing on loading package lattice

[David Winsemius]

> On Jan 28, 2018, at 4:50 AM, Michael Lane  wrote:
> 
> Hi R Help Team
> 
> I recently updated my R installation to R 3.4.3 and updated to later version 
> of R Studio and I found that the package semPLS will not load even though 
> installed and it seems to be failing on loading package lattice
> 
> Getting the following error message:
> 
> library(semPLS)
> Loading required package: lattice
> Error: package or namespace load failed for 'lattice':
> .onLoad failed in loadNamespace() for 'grid', details:
>  call: fun(libname, pkgname)
>  error: object 'C_initGrid' not found
> Error: package 'lattice' could not be loaded

Although it is the attempt to load lattice that triggers the fault, the 
difficulty actually reported was with loading the 'grid' Namespace. The grid 
package is one of the core R packages, so I have two suggestions:

First try to load grid as a package. (Although it's a core R package, it is not 
automatically loaded when R starts.)

library(grid)

If that fails then you clearly have a broken installation of R, and would need 
to reinstall. If it doesn't fail, then the error message would seem to be a 
puzzle, and you should re-read the Posting Guide and follow the instructions on 
reported unexpected behavior which you have so far only followed very 
incompletely.

-- 
David.


> 
> Any advice or help on this bug would be much appreciated
> 
> Best regards Michael
> 
> Dr Michael Lane USQ 
> Profile<http://staffprofile.usq.edu.au/Profile/Michael-Lane>
> PhD Information Systems, USQ
> Email: michael.l...@usq.edu.au<mailto:michael.l...@usq.edu.au>
> Ph 07 4631 1268
> Mobile 0407 316 391
> Academic Coordinator School of Management and Enterprise
> Member of Editoral Board Australasian Journal of Information 
> Systems<http://journal.acs.org.au/index.php/ajis>
> Member of Editoral Board Journal of Information Systems 
> Education<http://jise.org/>
> 
> 
> 
> _
> This email (including any attached files) is confidentia...{{dropped:18}}
> 
> __
> 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.

David Winsemius
Alameda, CA, USA

'Any technology distinguishable from magic is insufficiently advanced.'   
-Gehm's Corollary to Clarke's Third Law

__
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] Hausman test

[David Winsemius]

> On Feb 11, 2018, at 8:29 AM, PAOLO PILI  wrote:
> 
> you are right about the 3rd line but it doesn't help me for my problem. I
> remove the 3rd line but there is still the same problem:
> 
> Error in solve.default (dvcov):
>   the system is numerically unique: reciprocity condition value =
> 1.63418e-19

That suggests inclusion of too many categorical (factor) variables relative to 
the sample size in the predictor variables. Use tabular methods to investigate. 
Unable to be more specific in the absence of a proper description of the data 
situation.

-- 
David.
> 
> Paolo
> 
> 2018-02-11 16:54 GMT+01:00 Bert Gunter :
> 
>> Note the typo in your 3rd line: data <
>> 
>> Don't  know if this means anything...
>> 
>> Bert
>> 
>> 
>> 
>> On Feb 11, 2018 7:33 AM, "PAOLO PILI"  wrote:
>> 
>>> Hello,
>>> 
>>> I have a problem with Hausman test. I am performing my analysis with these
>>> commands:
>>> 
>>>> library(plm)
>>>> data<-read.csv2("paolo.csv",header=TRUE)
>>>> data<
>>> pdata.frame(data,index=c("FIRM","YEAR"),drop.index=TRUE,row.names=TRUE)
>>>> 
>>> RECEIV~LSIZE+LAGE+LAGE2+CFLOW+STLEV+FCOST+PGROWTH+NGROWTH+TU
>>> RN+GPROF+GPROF2
>>>> 
>>> grun.fe<-plm(RECEIV~LSIZE+LAGE+LAGE2+CFLOW+STLEV+FCOST+PGROW
>>> TH+NGROWTH+TURN+GPROF+GPROF2,data=data,model="within")
>>>> grun.re
>>> <-plm(RECEIV~LSIZE+LAGE+LAGE2+CFLOW+STLEV+FCOST+PGROWTH+NGRO
>>> WTH+TURN+GPROF+GPROF2,data=data,model="random")
>>>> 
>>> gw<-plm(RECEIV~LSIZE+LAGE+LAGE2+CFLOW+STLEV+FCOST+PGROWTH+
>>> NGROWTH+TURN+GPROF+GPROF2,data=data,model="within")
>>>> 
>>> gr<-plm(RECEIV~LSIZE+LAGE+LAGE2+CFLOW+STLEV+FCOST+PGROWTH+
>>> NGROWTH+TURN+GPROF+GPROF2,data=data,model="random")
>>>> phtest(gw,gr)
>>> 
>>> I got this answer:
>>> 
>>> Error in solve.default(dvcov) :
>>> 
>>> how can I solve this problem?
>>> 
>>> Thank you
>>> 
>>>[[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/posti
>>> ng-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.

David Winsemius
Alameda, CA, USA

'Any technology distinguishable from magic is insufficiently advanced.'   
-Gehm's Corollary to Clarke's Third Law

__
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] Parallel assignments and goto

[David Winsemius]

torial_tr_1(100),
>   factorial_tr_2(100),
>   factorial_tr_3(100))
> Unit: microseconds
>expr   min lqmean medianuq 
>max neval
>  factorial(100)52.47960.264093.4306967.513083.925   
> 2062.481   100
> factorial_tr_1(100) 8.875 9.652549.1959510.694511.217   
> 3818.823   100
> factorial_tr_2(100)  5296.350  5525.0745  5973.77664  5737.8730  6260.128   
> 8471.301   100
> factorial_tr_3(100) 77554.457 80757.0905 87307.28737 84004.0725 89859.169 
> 171039.228   100
> 
> I can live with the “introducing extra variables” solution to parallel 
> assignment, and I could hack my way out of using `with` or `bind` in 
> rewriting `cases`, but restarting a `repeat` loop would really make for a 
> nicer solution. I know that `goto` is considered harmful, but really, in this 
> case, it is what I want.
> 
> A `callCC` version also solves the problem
> 
> factorial_tr_4 <- function(n, acc = 1) {
>function_body <- function(continuation) {
>if (n <= 1) {
>continuation(acc)
>} else {
>continuation(list("continue", n = n - 1, acc = acc * n))
>}
>}
>repeat {
>result <- callCC(function_body)
>if (is.list(result) && result[[1]] == "continue") {
>n <- result$n
>acc <- result$acc
>next
>} else {
>return(result)
>}
>}
> }
> 
> But this requires that I know how to distinguish between a valid return value 
> and a tag for “next” and is still a lot slower than the `next` solution
> 
> microbenchmark::microbenchmark(factorial(100),
>   factorial_tr_1(100),
>   factorial_tr_2(100),
>   factorial_tr_3(100),
>   factorial_tr_4(100))
> Unit: microseconds
>expr   min lqmean medianuq 
>max neval
>  factorial(100)54.10961.809581.3316781.878589.748
> 243.554   100
> factorial_tr_1(100) 9.025 9.903511.3860711.199012.008 
> 22.375   100
> factorial_tr_2(100)  5272.524  5798.3965  6302.40467  6077.7180  6492.959   
> 9967.237   100
> factorial_tr_3(100) 66186.080 72336.2810 76480.75172 73632.9665 75405.054 
> 203785.673   100
> factorial_tr_4(100)   270.978   302.7890   337.48763   313.9930   334.096   
> 1425.702   100
> 
> I don’t necessarily need the tail-recursion optimisation to be faster than 
> the recursive version; just getting out of the problem of too deep recursions 
> is a benefit, but I would rather not pay with an order of magnitude for it. I 
> could, of course, try to handle cases that works with `next` in one way, and 
> other cases using `callCC`, but I feel it should be possible with a version 
> that handles all cases the same way.
> 
> Is there any way to achieve this?
> 
> Cheers
>   Thomas

I didn't see any reference to the R `Recall` or `local` functions. I don't 
remember that tail optimization is something that R provides, however.


David Winsemius
Alameda, CA, USA

'Any technology distinguishable from magic is insufficiently advanced.'   
-Gehm's Corollary to Clarke's Third Law

__
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] Help with regular expressions

[David Winsemius]

> On Feb 12, 2018, at 6:22 PM, Dennis Fisher  wrote:
> 
> R 3.4.2
> OS X
> 
> Colleagues
> 
> I would appreciate some help with regular expressions.
> 
> I have string that looks like:
>   " ITERATION  ,THETA1 ,THETA2
>  ,THETA3 ,THETA4 ,THETA5 
> ,THETA6 ,THETA7 ,SIGMA(1,1) 
> ,SIGMA(2,1) ,SIGMA(2,2)”
> 
> In the entries that contain:
>   (X,Y)   # for example, SIGMA(1,1)
> I would like to replace the comma with a period, e.g., SIGMA(1.1) but NOT the 
> other commas
> 
> The end-result would be:
>   " ITERATION  ,THETA1 ,THETA2
>  ,THETA3 ,THETA4 ,THETA5 
> ,THETA6 ,THETA7 ,SIGMA(1.1) 
> ,SIGMA(2.1) ,SIGMA(2.2)”
> 
> Can someone provide the regular expression code to accomplish this?

gsub( "([(]\\d+)([,])(\\d+[)])", "\\1.\\3", x)\
#---
[1] "ITERATION  ,THETA1 ,THETA2 
,THETA3 ,THETA4 ,THETA5 ,THETA6 
,THETA7 ,SIGMA(1.1) ,SIGMA(2.1) 
,SIGMA(2.2)"


> Thanks.
> 
> Dennis
> 
> Dennis Fisher MD
> P < (The "P Less Than" Company)
> Phone / Fax: 1-866-PLessThan (1-866-753-7784)
> www.PLessThan.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.

David Winsemius
Alameda, CA, USA

'Any technology distinguishable from magic is insufficiently advanced.'   
-Gehm's Corollary to Clarke's Third Law

__
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] Fleming-Harrington weighted log rank test

[David Winsemius]

> On Feb 13, 2018, at 4:02 PM, array chip via R-help  
> wrote:
> 
> Hi all, 
> 
> The survdiff() from survival package has an argument "rho" that implements 
> Fleming-Harrington weighted long rank test. 
> 
> But according to several sources including "survminer" package 
> (https://cran.r-project.org/web/packages/survminer/vignettes/Specifiying_weights_in_log-rank_comparisons.html),
>  Fleming-Harrington weighted log-rank test should have 2 parameters "p" and 
> "q" to control the weighting for earlier vs later times in the follow-up.
> 
> For example, setting rho=1 in survdiff() uses the Peto-Peto modification of 
> Gehan-Wilcox weights, which I can confirm by setting p=1 & 1=0 in comp() from 
> survminer package. similarly rho=0 is equivalent to p=0 & q=0
> 
> I am interested in putting more weights on survival difference in later 
> follow-up time. According to comp() from survminer package, that would set 
> p=0 & q=1 for Fleming-Harrington weights. 
> 
> My question is how I can do the same by setting certain values for "rho" in 
> the regular survival() function?

I think that survdiff uses a different version than what you have found. The 
G-rho family weights are:

w_j = [Sˆ(tj)]^ρ

So rather than two parameters on S(t) and (1-S(t)) as in the p,q version, you 
only have one parameter applied to S(t). This class handout says that the 
G-rho,gamma weighting scheme is not available in survdiff.

-- 
David Winsemius
Alameda, CA, USA

'Any technology distinguishable from magic is insufficiently advanced.'   
-Gehm's Corollary to Clarke's Third Law

__
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] Fleming-Harrington weighted log rank test

[David Winsemius]

> On Feb 14, 2018, at 5:26 PM, David Winsemius  wrote:
> 
>> 
>> On Feb 13, 2018, at 4:02 PM, array chip via R-help  
>> wrote:
>> 
>> Hi all, 
>> 
>> The survdiff() from survival package has an argument "rho" that implements 
>> Fleming-Harrington weighted long rank test. 
>> 
>> But according to several sources including "survminer" package 
>> (https://cran.r-project.org/web/packages/survminer/vignettes/Specifiying_weights_in_log-rank_comparisons.html),
>>  Fleming-Harrington weighted log-rank test should have 2 parameters "p" and 
>> "q" to control the weighting for earlier vs later times in the follow-up.
>> 
>> For example, setting rho=1 in survdiff() uses the Peto-Peto modification of 
>> Gehan-Wilcox weights, which I can confirm by setting p=1 & 1=0 in comp() 
>> from survminer package. similarly rho=0 is equivalent to p=0 & q=0
>> 
>> I am interested in putting more weights on survival difference in later 
>> follow-up time. According to comp() from survminer package, that would set 
>> p=0 & q=1 for Fleming-Harrington weights. 
>> 
>> My question is how I can do the same by setting certain values for "rho" in 
>> the regular survival() function?
> 
> I think that survdiff uses a different version than what you have found. The 
> G-rho family weights are:
> 
> w_j = [Sˆ(tj)]^ρ
> 
> So rather than two parameters on S(t) and (1-S(t)) as in the p,q version, you 
> only have one parameter applied to S(t). This class handout says that the 
> G-rho,gamma weighting scheme is not available in survdiff.
> 

Forgot to paste the link:

http://www.ics.uci.edu/~dgillen/STAT255/Handouts/lecture4.pdf

> -- 
> David Winsemius
> Alameda, CA, USA
> 
> 'Any technology distinguishable from magic is insufficiently advanced.'   
> -Gehm's Corollary to Clarke's Third Law

David Winsemius
Alameda, CA, USA

'Any technology distinguishable from magic is insufficiently advanced.'   
-Gehm's Corollary to Clarke's Third Law

__
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] RExcel issues

[David Winsemius]

> On Mar 1, 2018, at 2:02 PM, Michael Ashton  
> wrote:
> 
> Hi -
> 
> For a while I've used RExcel without problems to run a repeating portfolio 
> optimization problem where I solve for a portfolio allocation targeting a 
> particular risk, then solve for a different risk, etc. I call the commands 
> with (e.g.) rinterface.Rrun "(R command)"
> 
> Recently that macro started blowing up, returning #RErrors, and when I try to 
> trace the error I find that it is some kind of OLE error. I can't seem to 
> find anything about it online and it's hard to replicate for someone else 
> since you'd have to have my installation and my spreadsheet.
> 
> But I thought I'd ask the community and see if anyone else has had this 
> problem recently with Rexcel. I'm calling solnp from within the RRun command, 
> but again...it worked for a long time and I don't think it's the inputs that 
> are wonky. Perhaps an update of R was incompatible with Rexcel? I've updated 
> Rexcel to 3.2.16 but it made no difference.

I'm pretty sure that RExcel is a commercial program and I assume that licensed 
users are expected to bring problems to the vendor.

http://rcom.univie.ac.at/contact.html


> 
> Any suggestions of what to try will be warmly entertained!
> 
> Thanks,
> 
> Mike
> 
> Michael Ashton, CFA
> Managing Principal
> 
> Enduring Investments LLC
> W: 973.457.4602
> C: 551.655.8006

I assume with that address you have a commercial license.

-- 

David Winsemius
Alameda, CA, USA

'Any technology distinguishable from magic is insufficiently advanced.'   
-Gehm's Corollary to Clarke's Third Law

__
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] lmrob gives NA coefficients

[David Winsemius]

> On Mar 3, 2018, at 3:04 PM, Christien Kerbert  
> wrote:
> 
> Dear list members,
> 
> I want to perform an MM-regression. This seems an easy task using the
> function lmrob(), however, this function provides me with NA coefficients.
> My data generating process is as follows:
> 
> rho <- 0.15  # low interdependency
> Sigma <- matrix(rho, d, d); diag(Sigma) <- 1
> x.clean <- mvrnorm(n, rep(0,d), Sigma)

Which package are you using for mvrnorm?

> beta <- c(1.0, 2.0, 3.0, 4.0)
> error <- rnorm(n = n, mean = 0, sd = 1)
> y <- as.data.frame(beta[1]*rep(1, n) + beta[2]*x.clean[,1] +
> beta[3]*x.clean[,2] + beta[4]*x.clean[,3] + error)
> xy.clean <- cbind(x.clean, y)
> colnames(xy.clean) <- c("x1", "x2", "x3", "y")
> 
> Then, I pass the following formula to lmrob: f <- y ~ x1 + x2 + x3
> 
> Finally, I run lmrob: lmrob(f, data = data, cov = ".vcov.w")
> and this results in NA coefficients.

It would also be more courteous to specify the package where you are getting 
lmrob.

> 
> It would be great if anyone can help me out. Thanks in advance.
> 
> Regards,
> Christien
> 
>   [[alternative HTML version deleted]]

This is a plain text mailing list although it doesn't seem to have created 
problems this time.

> 
> __
> 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.

David Winsemius
Alameda, CA, USA

'Any technology distinguishable from magic is insufficiently advanced.'   
-Gehm's Corollary to Clarke's Third Law

__
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] Random effect in GAM (mgcv)

[David Winsemius]

> On Mar 4, 2018, at 7:37 AM, Bert Gunter  wrote:
> 
> Statistics questions are largely off  topic on this list, although they do
> sometimes intersect R programming issues, which are on topic. However, I
> believe a statistics list like stats.stackexchange.com might be more
> suitable for your query.
> 
> Cheers,
> Bert

What Bert says is certainly true of rhelp, but is not the case for the R mixed 
models mailing list. I think you should repost your question there. You should 
first read their mailing list info and configure you mail client to send plain 
text. Your current posting is marred by the presence of extraneous symbols that 
are the plain-text translation of HTML formatting. Those doubles asterisks are 
surely not in your formula. Furthermore, I for one do not see how one could 
determine "correctness" of an answer without access to the data.

Best of luck;
David.
> 
> 
> 
> Bert Gunter
> 
> "The trouble with having an open mind is that people keep coming along and
> sticking things into it."
> -- Opus (aka Berkeley Breathed in his "Bloom County" comic strip )
> 
> On Sun, Mar 4, 2018 at 4:10 AM, Vaniscotte Amélie 
> wrote:
> 
>> Dear R users,
>> 
>> 
>> I am using the *mgcv package* to model the ratio of hectares of damaged
>> culture by wild boar in french departments according to some
>> environmental covariates. I used a _Beta distribution_ for the response.
>> 
>> 
>> For each department, we estimated the damaged in 3 different culture
>> types (« Culture »). Our statistical individual are therefore the
>> department crossed by the culture type.
>> 
>> Also, the department are clustered into landscape types (« Cluster »).
>> 
>> Since I want to get the effect of the Culture type and the Landscape, I
>> keep those variables as fixed effects in the model.
>> 
>> 
>> Also, since we have 5 repetitions of the response and of some covariates
>> measurement in time per department and culture type, I put a random
>> effect on the Department per Culture type and the Year as fixed effect
>> as well.
>> 
>> 
>> The model takes the form :
>> 
>> 
>> *gam_tot <- gam (resp ~ Culture + Clust**er**:Culture + s(**Year**,k=4,
>> by=Culture) + s(**X1**, by=Culture) + s(**X2**, by=Culture) + s(Depts,
>> bs="re", by=Culture) *
>> 
>> *, family=betar(link="logit"),method="REML",data=data,select=FALSE)*
>> 
>> 
>> Then, I estimated the part of the model explained deviance provided by
>> each covariate. For that, I run the model without the given covariate
>> (keeping smooth parameters constant between models), and compute the
>> difference in deviance between the Full model (with the given covariate)
>> and the penalized model (without the given covariate):
>> 
>> (Full model Deviance – Penalized model Deviance) / Full Model Deviance
>> 
>> 
>> From that, I get a _huge proportion of Deviance explained by the random
>> effect_ (Department) of about 30 %, while the others covariates
>> explained less than 1 %.
>> 
>> 
>> 
>> *At this point, I have few questions :*
>> 
>> 
>> *- Do you think my model formula is correct regarding my data and
>> questions ?*
>> 
>> 
>> *- Is my estimate of explained deviance correct ?*
>> 
>> *In that case, how can I explain such a huge discrepancy between **the
>> part of deviance explained by **random and fixed effects ? *
>> 
>> 
>> Thanks for your help,
>> 
>> 
>> 
>> Amélie
>> 
>> 
>>[[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.

David Winsemius
Alameda, CA, USA

'Any technology distinguishable from magic is insufficiently advanced.'   
-Gehm's Corollary to Clarke's Third Law

__
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] data analysis for partial two-by-two factorial design

[David Winsemius]

> On Mar 5, 2018, at 8:52 AM, Ding, Yuan Chun  wrote:
> 
> Hi Bert,
> 
> I am very sorry to bother you again.
> 
> For the following question, as you suggested, I posted it in both Biostars 
> website and stackexchange website, so far no reply.
> 
> I really hope that you can do me a great favor to share your points about how 
> to explain the coefficients for drug A and drug B if run anova model 
> (response variable = drug A + drug B). is it different from running three 
> separate T tests?
> 
> Thank you so much!!
> 
> Ding
> 
> I need to analyze data generated from a partial two-by-two factorial design: 
> two levels for drug A (yes, no), two levels for drug B (yes, no);  however, 
> data points are available only for three groups, no drugA/no drugB, yes 
> drugA/no drugB, yes drugA/yes drug B, omitting the fourth group of no 
> drugA/yes drugB.  I think we can not investigate interaction between drug A 
> and drug B, can I still run  model using R as usual:  response variable = 
> drug A + drug B?  any suggestion is appreciated.

Replied on CrossValidated where this would be on-topic.

-- 
David,

> 
> 
> From: Bert Gunter [mailto:bgunter.4...@gmail.com]
> Sent: Friday, March 02, 2018 12:32 PM
> To: Ding, Yuan Chun
> Cc: r-help@r-project.org
> Subject: Re: [R] data analysis for partial two-by-two factorial design
> 
> 
> [Attention: This email came from an external source. Do not open attachments 
> or click on links from unknown senders or unexpected emails.]
> 
> 
> This list provides help on R programming (see the posting guide linked below 
> for details on what is/is not considered on topic), and generally avoids 
> discussion of purely statistical issues, which is what your query appears to 
> be. The simple answer is yes, you can fit the model as described,  but you 
> clearly need the off topic discussion as to what it does or does not mean. 
> For that, you might try the 
> stats.stackexchange.com statistical site.
> 
> Cheers,
> Bert
> 
> 
> Bert Gunter
> 
> "The trouble with having an open mind is that people keep coming along and 
> sticking things into it."
> -- Opus (aka Berkeley Breathed in his "Bloom County" comic strip )
> 
> On Fri, Mar 2, 2018 at 10:34 AM, Ding, Yuan Chun 
> mailto:ycd...@coh.org>> wrote:
> Dear R users,
> 
> I need to analyze data generated from a partial two-by-two factorial design: 
> two levels for drug A (yes, no), two levels for drug B (yes, no);  however, 
> data points are available only for three groups, no drugA/no drugB, yes 
> drugA/no drugB, yes drugA/yes drug B, omitting the fourth group of no 
> drugA/yes drugB.  I think we can not investigate interaction between drug A 
> and drug B, can I still run  model using R as usual:  response variable = 
> drug A + drug B?  any suggestion is appreciated.
> 
> Thank you very much!
> 
> Yuan Chun Ding
> 
> 
> -
> -SECURITY/CONFIDENTIALITY WARNING-
> This message (and any attachments) are intended solely...{{dropped:31}}

__
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] data analysis for partial two-by-two factorial design

[David Winsemius]

> On Mar 5, 2018, at 2:27 PM, Bert Gunter  wrote:
> 
> David:
> 
> I believe your response on SO is incorrect. This is a standard OFAT (one 
> factor at a time) design, so that assuming additivity (no interactions), the 
> effects of drugA and drugB can be determined via the model you rejected:

>> three groups, no drugA/no drugB, yes drugA/no drugB, yes drugA/yes drug B, 
>> omitting the fourth group of no drugA/yes drugB.

> 
> For example, if baseline control (no drugs) has a response of 0, drugA has an 
> effect of 1, drugB has an effect of 2, and the effects are additive, with no 
> noise we would have:
> 
> > d <- data.frame(drugA = c("n","y","y"),drugB = c("n","n","y"))

d2 <- data.frame(trt = c("Baseline","DrugA_only","DrugA_drugB")
> 
> > y <- c(0,1,3)
> 
> And a straighforward inear model recovers the effects:
> 
> > lm(y ~ drugA + drugB, data=d)
> 
> Call:
> lm(formula = y ~ drugA + drugB, data = d)
> 
> Coefficients:
> (Intercept)   drugAy   drugBy  
>   1.282e-161.000e+002.000e+00  

I think the labeling above is rather to mislead since what is labeled drugB is 
actually A&B. I think the method I suggest is more likely to be interpreted 
correctly:

> d2 <- data.frame(trt = c("Baseline","DrugA_only","DrugA_drugB"))
>  y <- c(0,1,3)
> lm(y ~ trt, data=d2)

Call:
lm(formula = y ~ trt, data = d2)

Coefficients:
   (Intercept)  trtDrugA_drugB   trtDrugA_only  
 2.564e-16   3.000e+00   1.000e+00  

-- 
David.
> 
> As usual, OFAT designs are blind to interactions, so that if they really 
> exist, the interpretation as additive effects is incorrect.
> 
> Cheers,
> Bert
> 
> 
> Bert Gunter
> 
> "The trouble with having an open mind is that people keep coming along and 
> sticking things into it."
> -- Opus (aka Berkeley Breathed in his "Bloom County" comic strip )
> 
> On Mon, Mar 5, 2018 at 2:03 PM, David Winsemius  
> wrote:
> 
> > On Mar 5, 2018, at 8:52 AM, Ding, Yuan Chun  wrote:
> >
> > Hi Bert,
> >
> > I am very sorry to bother you again.
> >
> > For the following question, as you suggested, I posted it in both Biostars 
> > website and stackexchange website, so far no reply.
> >
> > I really hope that you can do me a great favor to share your points about 
> > how to explain the coefficients for drug A and drug B if run anova model 
> > (response variable = drug A + drug B). is it different from running three 
> > separate T tests?
> >
> > Thank you so much!!
> >
> > Ding
> >
> > I need to analyze data generated from a partial two-by-two factorial 
> > design: two levels for drug A (yes, no), two levels for drug B (yes, no);  
> > however, data points are available only for three groups, no drugA/no 
> > drugB, yes drugA/no drugB, yes drugA/yes drug B, omitting the fourth group 
> > of no drugA/yes drugB.  I think we can not investigate interaction between 
> > drug A and drug B, can I still run  model using R as usual:  response 
> > variable = drug A + drug B?  any suggestion is appreciated.
> 
> Replied on CrossValidated where this would be on-topic.
> 
> --
> David,
> 
> >
> >
> > From: Bert Gunter [mailto:bgunter.4...@gmail.com]
> > Sent: Friday, March 02, 2018 12:32 PM
> > To: Ding, Yuan Chun
> > Cc: r-help@r-project.org
> > Subject: Re: [R] data analysis for partial two-by-two factorial design
> >
> > 
> > [Attention: This email came from an external source. Do not open 
> > attachments or click on links from unknown senders or unexpected emails.]
> > 
> >
> > This list provides help on R programming (see the posting guide linked 
> > below for details on what is/is not considered on topic), and generally 
> > avoids discussion of purely statistical issues, which is what your query 
> > appears to be. The simple answer is yes, you can fit the model as 
> > described,  but you clearly need the off topic discussion as to what it 
> > does or does not mean. For that, you might try the 
> > stats.stackexchange.com<http://stats.stackexchange.com> statistical site.
> >
> > Cheers,
> > Bert
> >
> >
> > Bert Gunter
> >
> > "The trouble with having an open mind is that people keep coming along and 
> > sticking things into it."
> > -- Opus (aka Berkeley Breathed in his "Bloom County" comic strip )
&g

Re: [R] data analysis for partial two-by-two factorial design

[David Winsemius]

> On Mar 5, 2018, at 3:04 PM, Bert Gunter  wrote:
> 
> But of course the whole point of additivity is to decompose the combined 
> effect as the sum of individual effects.

Agreed. Furthermore your encoding of the treatment assignments has the 
advantage that the default treatment contrast for A+B will have a statistical 
estimate associated with it. That was a deficiency of my encoding that Ding 
found problematic. I did have the incorrect notion that the encoding of Drug B 
in the single drug situation would have been NA and that the `lm`-function 
would produce nothing useful. Your setup had not occurred to me.

Best;
David.

> 
> "Mislead" is a subjective judgment, so no comment. The explanation I provided 
> is standard. I used it for decades when I taught in industry.
> 
> Cheers,
> Bert
> 
> 
> 
> Bert Gunter
> 
> "The trouble with having an open mind is that people keep coming along and 
> sticking things into it."
> -- Opus (aka Berkeley Breathed in his "Bloom County" comic strip )
> 
> On Mon, Mar 5, 2018 at 3:00 PM, David Winsemius  
> wrote:
> 
> > On Mar 5, 2018, at 2:27 PM, Bert Gunter  wrote:
> >
> > David:
> >
> > I believe your response on SO is incorrect. This is a standard OFAT (one 
> > factor at a time) design, so that assuming additivity (no interactions), 
> > the effects of drugA and drugB can be determined via the model you rejected:
> 
> >> three groups, no drugA/no drugB, yes drugA/no drugB, yes drugA/yes drug B, 
> >> omitting the fourth group of no drugA/yes drugB.
> 
> >
> > For example, if baseline control (no drugs) has a response of 0, drugA has 
> > an effect of 1, drugB has an effect of 2, and the effects are additive, 
> > with no noise we would have:
> >
> > > d <- data.frame(drugA = c("n","y","y"),drugB = c("n","n","y"))
> 
> d2 <- data.frame(trt = c("Baseline","DrugA_only","DrugA_drugB")
> >
> > > y <- c(0,1,3)
> >
> > And a straighforward inear model recovers the effects:
> >
> > > lm(y ~ drugA + drugB, data=d)
> >
> > Call:
> > lm(formula = y ~ drugA + drugB, data = d)
> >
> > Coefficients:
> > (Intercept)   drugAy   drugBy
> >   1.282e-161.000e+002.000e+00
> 
> I think the labeling above is rather to mislead since what is labeled drugB 
> is actually A&B. I think the method I suggest is more likely to be 
> interpreted correctly:
> 
> > d2 <- data.frame(trt = c("Baseline","DrugA_only","DrugA_drugB"))
> >  y <- c(0,1,3)
> > lm(y ~ trt, data=d2)
> 
> Call:
> lm(formula = y ~ trt, data = d2)
> 
> Coefficients:
>(Intercept)  trtDrugA_drugB   trtDrugA_only
>  2.564e-16   3.000e+00   1.000e+00
> 
> --
> David.
> >
> > As usual, OFAT designs are blind to interactions, so that if they really 
> > exist, the interpretation as additive effects is incorrect.
> >
> > Cheers,
> > Bert
> >
> >
> > Bert Gunter
> >
> > "The trouble with having an open mind is that people keep coming along and 
> > sticking things into it."
> > -- Opus (aka Berkeley Breathed in his "Bloom County" comic strip )
> >
> > On Mon, Mar 5, 2018 at 2:03 PM, David Winsemius  
> > wrote:
> >
> > > On Mar 5, 2018, at 8:52 AM, Ding, Yuan Chun  wrote:
> > >
> > > Hi Bert,
> > >
> > > I am very sorry to bother you again.
> > >
> > > For the following question, as you suggested, I posted it in both 
> > > Biostars website and stackexchange website, so far no reply.
> > >
> > > I really hope that you can do me a great favor to share your points about 
> > > how to explain the coefficients for drug A and drug B if run anova model 
> > > (response variable = drug A + drug B). is it different from running three 
> > > separate T tests?
> > >
> > > Thank you so much!!
> > >
> > > Ding
> > >
> > > I need to analyze data generated from a partial two-by-two factorial 
> > > design: two levels for drug A (yes, no), two levels for drug B (yes, no); 
> > >  however, data points are available only for three groups, no drugA/no 
> > > drugB, yes drugA/no drugB, yes drugA/yes drug B, omitting the fourth 
> > > group of no drugA/yes drugB.  I think we can not investigate interaction 
> > > between drug A and drug B, can I still run  model using R as usual:  
> > &g

Re: [R] raster time series statistics

[David Winsemius]

> On Mar 5, 2018, at 3:28 PM,  
>  wrote:
> 
> Hi List,
> 
> The following code returns an "Error in as.POSIXlt.character(x, tz, ...) :   
> character string is not in a standard unambiguous format"

I'm unable to produce that error. Which function was being evaluated to produce 
the error? I don't see where as.POSIXlt would have been called. You don't have 
any Date or POSIXt-classed variables. (I did need to also load the stringr 
package to get str_pad into my workspace.)

> 
> require(raster)
> require(rts)
> require(stringi)
> r <- raster(ncol=100, nrow=100)
> values(r) <- runif(ncell(r))
> list(ID=seq(1:24),month=rep(str_pad(1:12, pad = 0,width = 2 , 
> "left"),2),year=sort(rep(2016:2017,12)))->dt
> stack(r)->s
> r->rs
> for(i in 1:23){
> rs[]<-r[]*i
>  addLayer(s,rs)->s
> print(nlayers(s))
> }
> timelst<-paste0(unlist(dt['year']),'-',unlist(dt['month']))
> rts(s,time=as.yearmon(timelst))->rsts
> str(rsts@time)
> apply.monthly(rsts,mean)
> 
> I was expecting that the statistics accept the yearmonth format of the 
> timeslot, as it allows subsetting.
> Any suggestions?
> Thanks
> Herry
> 
>   [[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.

David Winsemius
Alameda, CA, USA

'Any technology distinguishable from magic is insufficiently advanced.'   
-Gehm's Corollary to Clarke's Third Law

__
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] raster time series statistics

[David Winsemius]

> On Mar 5, 2018, at 9:07 PM,  
>  wrote:
> 
> It works if you use as.Date. But this defeates the purpose for the yearmon 
> notion...
> 
> require(raster)
> require(rts)
> require(stringr)
> r <- raster(ncol=100, nrow=100)
> values(r) <- runif(ncell(r))
> stack(r)->s
> r->rs
> for(i in 1:23){
> rs[]<-r[]*i
>  addLayer(s,rs)->s
> print(nlayers(s))
> }
> dt<-list(ID=seq(1:24),month=rep(formatC(1:12,flag=0,width=2),2), 
> year=sort(rep(2016:2017,12)))
> timelst<-paste0(unlist(dt['year']),'-',unlist(dt['month']),"-01")
> strptime(timelst,format="%Y-%m-%d")->t1
> 
> 
> rts(s,time=as.yearmon(t1))->rsts
> subset(rsts,'2017')->r2017
> class(r2017@time)
> class(rsts@time)
> apply.monthly(rsts,mean) # this creates error
> 
> rts(s,time=as.Date(t1))->rsts1
> apply.monthly(rsts1,mean) # this creates output

So it appears that you have shown that apply.monthly doesn't work with 
'yearmon'-classed vectors in the endpoints of a RasterStackTS. It's not 
documented as doing so after all. So this is really a feature request for the 
maintainer. You should contact the responsible individual.

maintainer('rts')
[1] "Babak Naimi "

> 
> 
> -Original Message-
> From: Jim Lemon [mailto:drjimle...@gmail.com] 
> Sent: Tuesday, 6 March 2018 11:40 AM
> To: Herr, Alexander (L&W, Black Mountain) 
> Cc: David Winsemius ; r-help mailing list 
> 
> Subject: Re: [R] raster time series statistics
> 

David Winsemius
Alameda, CA, 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] subsetting comparison problem

[David Winsemius]

> On Mar 11, 2018, at 3:32 PM, Neha Aggarwal  wrote:
> 
> Hello All,
> I am facing a unique problem and am unable to find any help in R help pages
> or online. I will appreciate your help for the following problem:
> I have 2 data-frames, samples below and there is an expected output
> 
> R Dataframe1:
>C1  C2   C3 C4.. CN
> R1   0  1   0   1
> R21  0  11
> R31  0   0 0
> .
> .
> .
> RN
> 
> U Dataframe2 :
> C1 C2C3 C4.. CN
> U1 1   101
> U2 1   1 11
> 
> 
> Expected Output:
> U1 satisfies R1, R3
> U2 satisfies R1, R2, R3
> 

I don't think you have communicated what sort of meaning is attached to the 
word "satisfies".

Here's a double loop that reports membership of the column names of each row of 
U (Dataframe2) in each row of R (Dataframe1):

 apply( Dataframe2, 1, function(x){ z <- which(x==1);
   z2 <- names(x)[z];  
zlist=apply(Dataframe1, 1, function(y){ z3 <- 
which(y==1); 
z4 <- 
names(y)[z3]; 
z4[ which(z4 
%in% z2) ]}); 
zlist})
$U1
$U1$R1
[1] "C2" "C4"

$U1$R2
[1] "C1" "C4"

$U1$R3
[1] "C1"


$U2
$U2$R1
[1] "C2" "C4"

$U2$R2
[1] "C1" "C3" "C4"

$U2$R3
[1] "C1"

-- 
David.


> So this is a comparison of dataframes problem, with a subset dimension.
> There are 2 dataframe R and U. column names are same. There are certain
> columns belonging to each row in dataframe 1, denoted as 1s, while there
> are certain cols to each U denoted as 1s in each URow in dataframe2.
> 
> I have to find relationships between Rs and Us. So i start with each U row
> in U dataframe (lets say U1 row) and try to find all the rows in R
> dataframe, which are subset of U1 row.
> 
> I cant find a way to compare rows to see if one is subset of
> anotherwhat can I try, any pointers/ packages will be great help.
> Please help.
> 
> Thanks
> Neha
> 
>   [[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.

David Winsemius
Alameda, CA, USA

'Any technology distinguishable from magic is insufficiently advanced.'   
-Gehm's Corollary to Clarke's Third Law

__
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] Need Help on Batch forecasting for multiple items at one go

[David Winsemius]

> On Mar 12, 2018, at 8:47 PM, Manish Mukherjee  
> wrote:
> 
> Hi All,
> 
> 
> I have time series weekly data for number of servers for 62 weeks. i have to 
> forecast for next 8 weeks on what will be the usage.  i have a challenge in 
> the code where it is giving output for the last week value of  all the 
> servers, instead i need the output for next 8 weeks . i am attaching the data 
> and my r script for your reference.

Your attachment failed to come through. Your mail client failed to label it as 
an acceptable format.
> 
> 
> 
> Thanks & Regards
> Manish Mukherjee
> __
> 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.

David Winsemius
Alameda, CA, USA

'Any technology distinguishable from magic is insufficiently advanced.'   
-Gehm's Corollary to Clarke's Third Law

__
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] R project global options file

[David Winsemius]

> On Mar 15, 2018, at 10:24 PM, Adrian Friskin  
> wrote:
> 
> Hello R-help,
> I currently have R-project 3.4.2 and R-studio 1.1.383 
> installed on some of our universities computer labs. Since the installation 
> of Visual studio the default version of R installed has changed to the 
> version installed by VS. Users can change the default R version: found in R 
> studio global options but users need to do this every session. Is there a 
> file or setting found on the pc where this information is kept? If so, where 
> is it found.

This appears to be a question for the vendor of VS (guessing this would be 
Microsoft).
> 
> Cheers,
> 
> Adrian
> 
> Adrian Friskin | Senior Technology Support Officer | Learning Environments - 
> Labs Team | Learning Environments
> and Technology Services | Queensland University of Technology | Synergy 
> Building Level 3 - LETS, 88 Musk Avenue, Kelvin Grove 4059 |
> ph +61 (07) 3138 3941 | 
> adrian.fris...@qut.edu.au<mailto:adrian.fris...@qut.edu.au> | CRICOS No 00213J
> 
> 
>   [[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.

David Winsemius
Alameda, CA, USA

'Any technology distinguishable from magic is insufficiently advanced.'   
-Gehm's Corollary to Clarke's Third Law

__
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] Package 'pd.mirna.1.0.2xgain' was not found in the BioConductor repository

[David Winsemius]

> On Mar 20, 2018, at 11:49 PM, Persian Irani  wrote:
> 
> Hi all!
> 
> While I am trying to read .cel files with oligo package:
> 
> afbatch=read.celfiles(list.celfiles())
> 
> I get an error:
> Package 'pd.mirna.1.0.2xgain' was not found in the BioConductor repository
> 
> How can I overcome this?

I'm guessing that there is an error in the input file. That doesn't look like a 
valid package name. The BioC repo has a package names pd.mirna.1.0

Perhaps a comma or some other separator got misplaced along the way. At any 
rate the right place to go would be the BioC support pages, not Rhelp. I'm sure 
you will need to supply them with more information than you supplied here.

> 
> Thank you in advance
> 
> __
> 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.

David Winsemius
Alameda, CA, USA

'Any technology distinguishable from magic is insufficiently advanced.'   
-Gehm's Corollary to Clarke's Third Law

__
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] Sum of columns of a data frame equal to NA when all the elements are NA

[David Winsemius]

> On Mar 21, 2018, at 9:58 AM, Jeff Newmiller  wrote:
> 
> What do you mean by "should not"?
> 
> NULL means "missing object" in R. The result of the sum function is always 
> expected to be numeric... so NA_real or NA_integer could make sense as 
> possible return values. But you cannot compute on NULL so no, that doesn't 
> work. 
> 
> See the note under the "Value" section of ?sum as to why zero is returned 
> when all inputs are removed.

But your perspective above that sentence is at odds with the note in the 
Details section of the same document:

> sum(NULL)
[1] 0
> ?sum
starting httpd help server ... done
> sum( integer(0) )
[1] 0

Best;
David

> -- 
> Sent from my phone. Please excuse my brevity.
> 
> On March 21, 2018 9:03:29 AM PDT, Boris Steipe  
> wrote:
>> Should not the result be NULL if you have removed the NA with
>> na.rm=TRUE ?
>> 
>> B.
>> 
>> 
>> 
>>> On Mar 21, 2018, at 11:44 AM, Stefano Sofia
>>  wrote:
>>> 
>>> Dear list users,
>>> let me ask you this trivial question. I worked on that for a long
>> time, by now.
>>> Suppose to have a data frame with NAs and to sum some columns with
>> rowSums:
>>> 
>>> df <- data.frame(A = runif(10), B = runif(10), C = rnorm(10))
>>> df[1, ] <- NA
>>> rowSums(df[ , which(names(df) %in% c("A","B"))], na.rm=T)
>>> 
>>> If all the elements of the selected columns are NA, rowSums returns 0
>> while I need NA.
>>> Is there an easy and efficient way to use rowSums within a function
>> like
>>> 
>>> function(x) ifelse(all(is.na(x)), as.numeric(NA), rowSums...)?
>>> 
>>> or an equivalent function?
>>> 
>>> Thank you for your help
>>> Stefano
>>> 
>>> 
>>> 
>>>(oo)
>>> --oOO--( )--OOo
>>> Stefano Sofia PhD
>>> Area Meteorologica e  Area nivologica - Centro Funzionale
>>> Servizio Protezione Civile - Regione Marche
>>> Via del Colle Ameno 5
>>> 60126 Torrette di Ancona, Ancona
>>> Uff: 071 806 7743
>>> E-mail: stefano.so...@regione.marche.it
>>> ---Oo-oO
>>> 
>>> 
>>> 
>>> AVVISO IMPORTANTE: Questo messaggio di posta elettronica può
>> contenere informazioni confidenziali, pertanto è destinato solo a
>> persone autorizzate alla ricezione. I messaggi di posta elettronica per
>> i client di Regione Marche possono contenere informazioni confidenziali
>> e con privilegi legali. Se non si è il destinatario specificato, non
>> leggere, copiare, inoltrare o archiviare questo messaggio. Se si è
>> ricevuto questo messaggio per errore, inoltrarlo al mittente ed
>> eliminarlo completamente dal sistema del proprio computer. Ai sensi
>> dell’art. 6 della DGR n. 1394/2008 si segnala che, in caso di necessità
>> ed urgenza, la risposta al presente messaggio di posta elettronica può
>> essere visionata da persone estranee al destinatario.
>>> IMPORTANT NOTICE: This e-mail message is intended to be received only
>> by persons entitled to receive the confidential information it may
>> contain. E-mail messages to clients of Regione Marche may contain
>> information that is confidential and legally privileged. Please do not
>> read, copy, forward, or store this message unless you are an intended
>> recipient of it. If you have received this message in error, please
>> forward it to the sender and delete it completely from your computer
>> system.
>>> 
>>> --
>>> Questo messaggio  stato analizzato da Libra ESVA ed  risultato non
>> infetto.
>>> This message was scanned by Libra ESVA and is believed to be clean.
>>> 
>>> 
>>> [[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.

David Winsemius
Alameda, CA, USA

'Any technology distinguishable from magic is insufficiently advanced.'   
-Gehm's Corollary to Clarke's Third Law

__
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] MCMCglmm multinomial model results

[David Winsemius]

> On Mar 22, 2018, at 1:31 PM, Michelle Kline  
> wrote:
> 
> Hi,
> 
> Thanks in advance for any help on this question. I'm running multinomial
> models using the MCMCglmm package. The models have 5 outcome variables
> (each with count data), and an additional two random effects built into the
> models. The issue is that when I use the following code, the summary only
> gives me results for four of the outcome variables.
> 
> Here is the code for my model:
> 
> m3.random <- MCMCglmm(cbind(Opp_teacher , Dir_teacher, Enh_teacher,
> SocTol_teacher, Eval_teacher) ~ trait -1,
>   random = ~ us(trait):other + us(trait):focal,
>   rcov = ~ us(trait):units,
>   prior = list(
> R = list(fix=1, V=0.5 * (I + J), n = 4),
> G = list(
>   G1 = list(V = diag(4), n = 4),
>   G2 = list(V = diag(4), n = 4))),
>   burnin = burn,
>   nitt = iter,
>   family = "multinomial5",
>   data = data,

We have no way to debug this without the data. Perhaps you should contact the 
maintainer and in your message attach the data?

 maintainer('MCMCglmm')
[1] "Jarrod Hadfield "


An equally effective approach would be to post (again with data that reproduces 
the error)  on the R-SIG-mixed-models mailing list since Hadfield is a regular 
contributor on that list. (To me it suggests not an error since you got output 
but rather a warning. Generally warnings and errors are properly labeled so you 
may not have included the full output.) 

-- 
David.
>   pr=TRUE,
>   pl=TRUE,
>   DIC = TRUE,
>   verbose = FALSE)
> 
> And the summary of the main effects:
> 
> post.mean  l-95% CI  u-95% CI eff.samppMCMC
> traitOpp_teacher-3.828752 -4.616731 -3.067424 184.4305 5.263158e-05
> traitDir_teacher-3.400481 -4.041069 -2.813063 259.1084 5.263158e-05
> traitEnh_teacher-1.779129 -2.197415 -1.366496 624.9759 5.263158e-05
> traitSocTol_teacher -2.852684 -3.429799 -2.332909 468.7098 5.263158e-05
> 
> 
> It is not an issue of the suppressing the intercept, since I'm already
> doing that (see the -1 term. When I remove that term, the model solutions
> includes an intercept and only 3 additional main effects).
> 
> The model does throw the following error, but after searching previous
> messages on this list, I've concluded that this error message doesn't have
> to do with  my current problem. Just in case: " observations with zero
> weight not used for calculating dispersion"
> 
> I have also posted a similar question on stackoverflow about a week ago,
> but with no response, so I thought I would try here. Link in case people
> want to gain reputation points for a
> response: 
> https://stackoverflow.com/questions/49309027/missing-term-in-mcmcglmm-multinomial-model-results-not-in-intercept-issue
> <https://stackoverflow.com/questions/49309027/missing-term-in-mcmcglmm-multinomial-model-results-not-in-intercept-issue>
> 
> And of course I've checked various other sources including the course
> notes, but can't make sense of why the 5th term is dropped from the model.
> Any help is much appreciated.
> 
> Best,
> 
> Michelle
> 
> -- 
> Michelle A. Kline, PhD
> 
> Assistant Professor
> Department of Psychology
> Simon Fraser 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.

David Winsemius
Alameda, CA, USA

'Any technology distinguishable from magic is insufficiently advanced.'   
-Gehm's Corollary to Clarke's Third Law

__
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] Calculate weighted proportions for several factors at once

[David Winsemius]

> On Mar 22, 2018, at 3:34 PM, Striessnig, Erich  
> wrote:
> 
> Hi,
> 
> I have a grouped data set and would like to calculate weighted proportions 
> for a large number of factor variables within each group member. Rather than 
> using dplyr::count() on each of these factors individually, the idea would be 
> to do it for all factors at once. Does anyone know how this would work? Here 
> is a reproducible example:
> 
> 
> # reproducible example
> df1 <- data.frame(wt=rnorm(90),
>  group=paste0('reg', 1:5),
>  var1=rep(c('male','female'), times=45),
>  var2=rep(c('low','med','high'), each=30)) %>% tbl_df()
> 
> # instead of doing this separately for each factor ...
> df2 <- df1 %>%
>  group_by(group) %>%
>  dplyr::count(var1, wt=wt) %>%
>  mutate(prop1=n/sum(n))
> 
> df3 <- df1 %>%
>  group_by(group) %>%
>  dplyr::count(var2, wt=wt) %>%
>  mutate(prop2=n/sum(n)) %>%
>  left_join(df2, by='group')
> 
> # I would like to do something like the following (which does of course not 
> work):
> my_fun <- function(x,wt){
>  freq1 <- dplyr::count(x, wt=wt)
>  prop1 <- freq1 / sum(freq1)
>  return(prop)
> }
> 
> df1 %>%
>  group_by(group) %>%
>  summarise_all(.funs=my_fun(.), .vars=c('var1', 'var2'))
> 

You might find useful functions in the ‘freqweights’ package. It appears from 
its description that it was design to fit into the tidyverse paradigm. I think 
the survey package might also be useful, but it is not particularly designed 
for use with tibbles and `%>%`. Might work. Might not.

HTH;
Dadid.
 
> 
> Best regards,
> Erich
> 
>   [[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.

David Winsemius
Alameda, CA, USA

'Any technology distinguishable from magic is insufficiently advanced.'   
-Gehm's Corollary to Clarke's Third Law

__
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] Medical risk package calculation RSI

[David Winsemius]

> On Mar 25, 2018, at 7:53 AM, Ramesh YAPALPARVI  
> wrote:
> 
> Hi all,
> 
> 
> 
> I'm using the medical risk package to determine the risk stratification Index 
> based on the ICD9 codes. Although, I have been successful in using it, I'm 
> unable to interpret the output.
> 
> 
> 
> here is the sample code
> 
> 
> 
> # Calculate RSI for each patient ("id") in dataframe
> 
> 
> cases <- data.frame(id=c(1,1,1,2,2,2,2,2),
> + 
> icd9cm=c("D4019","D25000","DV707","D71945","DV4365","D78079","D70909","D1958"))
> 
> 
> library(plyr)
> ddply(cases, .(id), function(x) { icd9cm_sessler_rsi(x$icd9cm) } )
> 
> 
> 
> 
> Output:
> 
> id rsi_1yrpod rsi_30dlos rsi_30dpod rsi_inhosp
> 1  1 -0.2860474  1.2481208 -0.5722244 -0.2612804
> 2  2  0.1417417  0.9779779 -1.0318395  0.000
> 
> Could anyone please help me?

That is just the code and error-free output from the help page for 
`icd9cm_sessler_rsi`.

Your posting does not actually conform to the rhelp mailing list standards 
since there is no illustrated problem with the use of the R language and rhelp 
is not designed to fill in gaps in your understanding of health services 
research methods. Have you read through and more importantly "worked through" 
the two vignettes that ship with that package? The readers of this list might 
be more motivated if you had presented evidence that you had made a good faith 
effort with the material that had been already made available to you.

vignette(pac='medicalrisk')

Vignettes in package ‘medicalrisk’:

medicalrisk  medicalrisk: Calculating risk and comorbidities
 from ICD-9-CM codes (source, html)
advanced medicalrisk: Correlating ICU days with mortality
 risk and comorbidities (source, html)

Once you have made an effort to "do your own homework", you might be in a 
position to construct a more substantial and more specific question for a venue 
where this might be more on-topic, perhaps CrossValidated.com 
<https://stats.stackexchange.com/>
-- 

David Winsemius
Alameda, CA, USA

'Any technology distinguishable from magic is insufficiently advanced.'   
-Gehm's Corollary to Clarke's Third Law

__
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] Names of variables needed in newdata for predict.glm

[David Winsemius]

> On Mar 31, 2018, at 8:48 AM, Bendix Carstensen  
> wrote:
> 
> all.vars works fine, EXCEPT, it give a bit too much.
> I only want the regression variables, but in the following example I also get 
> "k" the variable holding the chosen knots. Any machinery to find only "real" 
> regression variables?
> cheers, Bendix
> 
> library( splines )
> y <- rnorm(100)
> x <- rnorm(100)
> k <- -1:1
> ml <-  lm( y ~ bs(x,knots=k) )
> mg <- glm( y ~ bs(x,knots=k) )
> all.vars(ml$terms)
> all.vars(mg$terms)
> all.vars(mg$formula)

If you allowed a requirement that "real" regression variables have been passed 
in a data argument, then this might succeed:

> ml <-  lm( y ~ bs(x,knots=k), data=dat )
> all.vars(ml$terms)
[1] "y" "x" "k"
> all.vars(ml$formula)
character(0)
> all.vars(ml$terms)[ all.vars(ml$terms) %in% names(dat)]
[1] "y" "x"

-- 
David.
> 

> 
> Fra: Marc Girondot 
> Sendt: 8. marts 2018 06:26
> Til: Bendix Carstensen; r-help@r-project.org
> Emne: Re: [R] Names of variables needed in newdata for predict.glm
> 
> Hi,
> 
> Some try:
>> names(mi$xlevels)
> [1] "f"
>> all.vars(mi$formula)
> [1] "D" "x" "f" "Y"
>> names(mx$xlevels)
> [1] "f"
>> all.vars(mx$formula)
> [1] "D" "x" "f"
> 
> When offset is indicated out of the formula, it does not work...
> 
> Marc
> 
> Le 07/03/2018 à 06:20, Bendix Carstensen a écrit :
>> I would like to extract the names, modes [numeric/factor] and levels
>> of variables needed in a data frame supplied as newdata= argument to
>> predict.glm()
>> 
>> Here is a small example illustrating my troubles; what I want from
>> (both of) the glm objects is the vector c("x","f","Y") and an
>> indication that f is a factor:
>> 
>> library( splines )
>> dd <- data.frame( D = sample(0:1,200,rep=T),
>>   x = abs(rnorm(200)),
>>   f = factor(sample(letters[1:4],200,rep=T)),
>>   Y = runif(200,0.5,10) )
>> mx <- glm( D ~ ns(x,knots=1:2,Bo=c(0,5)) + f:I(x^2) , offset=log(Y) , 
>> family=poisson, data=dd)
>> mi <- glm( D ~ ns(x,knots=1:2,Bo=c(0,5)) + f:I(x^2) + offset(log(Y)), 
>> family=poisson, data=dd)
>> 
>> attr(mx$terms,"dataClasses")
>> attr(mi$terms,"dataClasses")
>> mi$xlevels
>> mx$xlevels
>> 
>> ...so far not quite there.
>> 
>> Regards,
>> 
>> Bendix Carstensen
>> 
>> Senior Statistician
>> Steno Diabetes Center
>> Clinical Epidemiology
>> Niels Steensens Vej 2-4
>> DK-2820 Gentofte, Denmark
>> b...@bxc.dk
>> bendix.carsten...@regionh.dk
>> http://BendixCarstensen.com
>> 
>> 
>> 
>> 
>> Denne e-mail indeholder fortrolig information. Hvis du ikke er den rette 
>> modtager af denne e-mail eller hvis du modtager den ved en fejltagelse, 
>> beder vi dig venligst informere afsender om fejlen ved at bruge 
>> svarfunktionen. Samtidig bedes du slette e-mailen med det samme uden at 
>> videresende eller kopiere den.
>> 
>> __
>> 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.
>> 
> 
> 
> 
> 
> 
> Denne e-mail indeholder fortrolig information. Hvis du ikke er den rette 
> modtager af denne e-mail eller hvis du modtager den ved en fejltagelse, beder 
> vi dig venligst informere afsender om fejlen ved at bruge svarfunktionen. 
> Samtidig bedes du slette e-mailen med det samme uden at videresende eller 
> kopiere den.
> 
> __
> 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.

David Winsemius
Alameda, CA, USA

'Any technology distinguishable from magic is insufficiently advanced.'   
-Gehm's Corollary to Clarke's Third Law

__
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] Obtain gradient at multiple values for exponetial decay model

[David Winsemius]

> On Apr 5, 2018, at 2:00 PM, g l  wrote:
> 
>> Sent: Thursday, April 05, 2018 at 4:40 PM
>> From: "Jeff Newmiller" 
>> 
>> the coef function.
>> 
> 
> For the benefit of other novices, used the following command to read the 
> documentation:
> 
> ?coef
> 
> Then tried and obtained:
> 
>> cvalue100<-coef(graphmodelp~100)
>> cvalue100
> NULL

Should have been:

  coef(graphmodelp)


> 
> Then looked at the model values which of course correspond to original 
> non-modelled values.
> 
> graphmodelp
>1 2 3 4 5 6 7 
> 8 
> 91.244636 69.457794 52.873083 40.248368 30.638107 23.322525 17.753714 
> 13.514590 
>91011 
> 10.287658  7.831233  5.961339

Read up on ?predict and what it delivers when only a model is offered as input.

> 
> This prompted to think that interpolation is required, but the function 
> 'approx' only seems to perform constant interpolation.
> 
> Is the correct thinking to find a function to perform interpolation, then 
> find/write a function to differentiate the model at a specific value of x, to 
> find gradient at that point?

Not correct. You already have `predict`. It is capale of using the `newdata` 
values to do interpolation with the values of the coefficients in the model. 
See: 

?predict

The original question asked for a derivative (i.e. a "gradient"), but so far 
it's not clear that you understand the mathematical definiton of that term. We 
also remain unclear whether this is homework.


> 
> __
> 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.

David Winsemius
Alameda, CA, USA

'Any technology distinguishable from magic is insufficiently advanced.'   
-Gehm's Corollary to Clarke's Third Law

__
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] Obtain gradient at multiple values for exponential decay model

[David Winsemius]

> On Apr 6, 2018, at 3:43 AM, g l  wrote:
> 
>> Sent: Friday, April 06, 2018 at 5:55 AM
>> From: "David Winsemius" 
>> 
>> 
>> Not correct. You already have `predict`. It is capale of using the `newdata` 
>> values to do interpolation with the values of the coefficients in the model. 
>> See: 
>> 
>> ?predict
>> 
> 
> The § details did not mention interpolation explicity; thanks.
> 
>> The original question asked for a derivative (i.e. a "gradient"), but so far 
>> it's not clear that you understand the mathematical definiton of that term. 
>> We also remain unclear whether this is homework.
>> 
> 
> The motivation of this post was simple differentiation of a tangent point 
> (dy/dx) manually, then wondering how to re-think in modern-day computing 
> terms. Hence the original question about asking the appropriate 
> functions/syntax to read further ("curiosity"), not the answer (indeed, 
> "homework"). :)
> 
> Personal curiosity should be considered "homework".

Besides symbolic differentiation, there is also the option of numeric 
differentiation. Here's an amateurish attempt:

myNumDeriv <- function(x){ (exp( predict (graphmodeld, 
newdata=data.frame(t=x+.0001))) - 
    exp( predict (graphmodeld, 
newdata=data.frame(t=x) )))/
  .0001 }
myNumDeriv(c(100, 250, 350))



David Winsemius
Alameda, CA, USA

'Any technology distinguishable from magic is insufficiently advanced.'   
-Gehm's Corollary to Clarke's Third Law

__
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.