[R] recode: how to avoid nested ifelse

[Paul Johnson]

In our Summer Stats Institute, I was asked a question that amounts to
reversing the effect of the contrasts function (reconstruct an ordinal
predictor from a set of binary columns). The best I could think of was to
link together several ifelse functions, and I don't think I want to do this
if the example became any more complicated.

I'm unable to remember a less error prone method :). But I expect you might.

Here's my working example code

## Paul Johnson pauljohn at ku.edu
## 2013-06-07

## We need to create an ordinal factor from these indicators
## completed elementary school
es - c(0, 0, 1, 0, 1, 0, 1, 1)
## completed high school
hs - c(0, 0, 1, 0, 1, 0, 1, 0)
## completed college graduate
cg - c(0, 0, 0, 0, 1, 0, 1, 0)

ed - ifelse(cg == 1, 3,
 ifelse(hs == 1, 2,
ifelse(es == 1, 1, 0)))

edf - factor(ed, levels = 0:3,  labels = c(none, es, hs, cg))
data.frame(es, hs, cg, ed, edf)

## Looks OK, but what if there are missings?
es - c(0, 0, 1, 0, 1, 0, 1, 1, NA, NA)
hs - c(0, 0, 1, 0, 1, 0, 1, 0, 1, NA)
cg - c(0, 0, 0, 0, 1, 0, 1, 0, NA, NA)
ed - ifelse(cg == 1, 3,
 ifelse(hs == 1, 2,
ifelse(es == 1, 1, 0)))
cbind(es, hs, cg, ed)

## That's bad, ifelse returns NA too frequently.
## Revise (becoming tedious!)

ed - ifelse(!is.na(cg)  cg == 1, 3,
 ifelse(!is.na(hs)  hs == 1, 2,
ifelse(!is.na(es)  es == 1, 1,
   ifelse(is.na(es), NA, 0
cbind(es, hs, cg, ed)


## Does the project director want us to worry about
## logical inconsistencies, such as es = 0 but cg = 1?
## I hope not.

Thanks in advance, I hope you are having a nice summer.

pj

-- 
Paul E. Johnson
Professor, Political Science  Assoc. Director
1541 Lilac Lane, Room 504  Center for Research Methods
University of Kansas University of Kansas
http://pj.freefaculty.org   http://quant.ku.edu

[[alternative HTML version deleted]]

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

Re: [R] recode: how to avoid nested ifelse

[Paul Johnson]

Thanks, guys.


On Sat, Jun 8, 2013 at 2:17 PM, Neal Fultz nfu...@gmail.com wrote:

 rowSums and Reduce will have the same problems with bad data you alluded
 to earlier, eg
 cg = 1, hs = 0

 But that's something to check for with crosstabs anyway.


This wrong data thing is a distraction here.  I guess I'd have to craft 2
solutions, depending on what the researcher says. (We can't assume es = 0
or es = NA and cg = 1 is bad data. There are some people who finish college
without doing elementary school (wasn't Albert Einstein one of those?) or
high school. I once went to an eye doctor who didn't finish high school,
but nonetheless was admitted to optometrist school.)

I did not know about the Reduce function before this. If we enforce the
ordering and clean up the data in the way you imagine, it would work.

I think the pmax is the most teachable and dependably not-getting-wrongable
approach if the data is not wrong.


 Side note: you should check out the microbenchmark pkg, it's quite handy.


Perhaps the working example of microbenchmark is the best thing in this
thread! I understand the idea behind it, but it seems like I can never get
it to work right. It helps to see how you do it.


 Rrequire(microbenchmark)
 Rmicrobenchmark(
 +   f1(cg,hs,es),
 +   f2(cg,hs,es),
 +   f3(cg,hs,es),
 +   f4(cg,hs,es)
 + )
 Unit: microseconds
expr   min lq median uq   max neval
  f1(cg, hs, es) 23029.848 25279.9660 27024.9640 29996.6810 55444.112   100
  f2(cg, hs, es)   730.665   755.5750   811.7445   934.3320  6179.798   100
  f3(cg, hs, es)85.029   101.6785   129.8605   196.2835  2820.187   100
  f4(cg, hs, es)   762.232   804.4850   843.7170  1079.0800 24869.548   100

 On Fri, Jun 07, 2013 at 08:03:26PM -0700, Joshua Wiley wrote:
  I still argue for na.rm=FALSE, but that is cute, also substantially
 faster
 
  f1 - function(x1, x2, x3) do.call(paste0, list(x1, x2, x3))
  f2 - function(x1, x2, x3) pmax(3*x3, 2*x2, es, 0, na.rm=FALSE)
  f3 - function(x1, x2, x3) Reduce(`+`, list(x1, x2, x3))
  f4 - function(x1, x2, x3) rowSums(cbind(x1, x2, x3))
 
  es - rep(c(0, 0, 1, 0, 1, 0, 1, 1, NA, NA), 1000)
  hs - rep(c(0, 0, 1, 0, 1, 0, 1, 0, 1, NA), 1000)
  cg - rep(c(0, 0, 0, 0, 1, 0, 1, 0, NA, NA), 1000)
 
  system.time(replicate(1000, f1(cg, hs, es)))
  system.time(replicate(1000, f2(cg, hs, es)))
  system.time(replicate(1000, f3(cg, hs, es)))
  system.time(replicate(1000, f4(cg, hs, es)))
 
   system.time(replicate(1000, f1(cg, hs, es)))
 user  system elapsed
22.730.03   22.76
   system.time(replicate(1000, f2(cg, hs, es)))
 user  system elapsed
 0.920.040.95
   system.time(replicate(1000, f3(cg, hs, es)))
 user  system elapsed
 0.190.020.20
system.time(replicate(1000, f4(cg, hs, es)))
 user  system elapsed
 0.950.030.98
 
 
  R version 3.0.0 (2013-04-03)
  Platform: x86_64-w64-mingw32/x64 (64-bit)
 
 
 
 
  On Fri, Jun 7, 2013 at 7:25 PM, Neal Fultz nfu...@gmail.com wrote:
   I would do this to get the highest non-missing level:
  
   x - pmax(3*cg, 2*hs, es, 0, na.rm=TRUE)
  
   rock chalk...
  
   -nfultz
  
   On Fri, Jun 07, 2013 at 06:24:50PM -0700, Joshua Wiley wrote:
   Hi Paul,
  
   Unless you have truly offended the data generating oracle*, the
   pattern: NA, 1, NA, should be a data entry error --- graduating HS
   implies graduating ES, no?  I would argue fringe cases like that
   should be corrected in the data, not through coding work arounds.
   Then you can just do:
  
   x - do.call(paste0, list(es, hs, cg))
  
table(factor(x, levels = c(000, 100, 110, 111), labels =
 c(none, es,hs, cg)))
   none   es   hs   cg
  4112
  
   Cheers,
  
   Josh
  
   *Drawn from comments by Judea Pearl one lively session.
  
  
   On Fri, Jun 7, 2013 at 6:13 PM, Paul Johnson pauljoh...@gmail.com
 wrote:
In our Summer Stats Institute, I was asked a question that amounts
 to
reversing the effect of the contrasts function (reconstruct an
 ordinal
predictor from a set of binary columns). The best I could think of
 was to
link together several ifelse functions, and I don't think I want to
 do this
if the example became any more complicated.
   
I'm unable to remember a less error prone method :). But I expect
 you might.
   
Here's my working example code
   
## Paul Johnson pauljohn at ku.edu
## 2013-06-07
   
## We need to create an ordinal factor from these indicators
## completed elementary school
es - c(0, 0, 1, 0, 1, 0, 1, 1)
## completed high school
hs - c(0, 0, 1, 0, 1, 0, 1, 0)
## completed college graduate
cg - c(0, 0, 0, 0, 1, 0, 1, 0)
   
ed - ifelse(cg == 1, 3,
 ifelse(hs == 1, 2,
ifelse(es == 1, 1, 0)))
   
edf - factor(ed, levels = 0:3,  labels = c(none, es, hs,
 cg))
data.frame(es, hs, cg, ed, edf)
   
## Looks OK, but what if there are missings?
es - c(0, 0, 1

[R] Windows R_LIBS_USER confusion under R-3.0.1

[Paul Johnson]

I would appreciate ideas about MS Windows install issues. I'm at our stats
summer camp and have been looking at a lot of Windows R installs and there
are some wrinkles about R_LIBS_USER.

On a clean Win7 or Win8 system, with R-3.0.1, we see the user library for
packages defaulting to  $HOME/R/win-library.

I think that's awesome, the way it should be. Yea! But it does not appear
that way on all systems, and I think it is because of lingering
after-effects of previous R installs.

In previous versions of R, R_LIBS_USER defaulted to
$HOME/AppData/Roaming/... That was not so great because it was (default)
hidden to the students and they were disoriented about how something that
does not exist (or show) could hold packages. Aside from teaching them how
to configure the file manager, we could navigate that.

The problem is that with R-3.0.1, sometimes we are seeing the user package
installs going into $HOME/AppData/Roaming/

In the user account, there is no $HOME/.Rprofile or $HOME/.Renviron.

I hate to tell non-expert users that they ought to go fishing in the
Windows registry, but I'm starting to suspect that is what they ought to
do.  What do you think?

PJ

-- 
Paul E. Johnson
Professor, Political Science  Assoc. Director
1541 Lilac Lane, Room 504  Center for Research Methods
University of Kansas University of Kansas
http://pj.freefaculty.org   http://quant.ku.edu

[[alternative HTML version deleted]]

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

Re: [R] 2SLS / TSLS / SEM non-linear

[Paul Johnson]

Please consider using the R package systemfit. It has existed for about 10
years, I've used it many times happily :)

I've not used gmm package, I'm not criticizing it. But I have good results
from systemfit. It has good documentation.

This package contains functions for fitting
 simultaneous systems of linear and nonlinear
 equations using Ordinary Least Squares (OLS),
 Weighted Least Squares (WLS), Seemingly Unrelated
 Regressions (SUR), Two-Stage Least Squares (2SLS),
 Weighted Two-Stage Least Squares (W2SLS), and
 Three-Stage Least Squares (3SLS).

If that doesn't work for you, please show us your example code, along with
the error messages.





On Sat, Jun 29, 2013 at 11:18 AM, John Fox j...@mcmaster.ca wrote:

 Dear Hans-Christian

  -Original Message-
  From: r-help-boun...@r-project.org [mailto:r-help-bounces@r-
  project.org] On Behalf Of hck
  Sent: Saturday, June 29, 2013 11:19 AM
  To: r-help@r-project.org
  Subject: Re: [R] 2SLS / TSLS / SEM non-linear
 
  The challenge is to firstly calculate the reduced form. As far as I
  know, the
  SEM package does not do this automatically. Am I correct?
 

 Right. If you look at the code in sem:::tsls.default, you'll see that it
 formulates and solves the 2SLS estimating equations directly (using a
 Cholesky decomposition). Moreover, tsls() is for linear structural
 equations.

 Best,
  John

 ---
 John Fox
 Senator McMaster Professor of Social Statistics
 Department of Sociology
 McMaster University
 Hamilton, Ontario, Canada


 
 
 
 
  --
  View this message in context: http://r.789695.n4.nabble.com/2SLS-TSLS-
  SEM-non-linear-tp4670123p4670595.html
  Sent from the R help mailing list archive at Nabble.com.
 
  __
  R-help@r-project.org mailing list
  https://stat.ethz.ch/mailman/listinfo/r-help
  PLEASE do read the posting guide http://www.R-project.org/posting-
  guide.html
  and provide commented, minimal, self-contained, reproducible code.

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




-- 
Paul E. Johnson
Professor, Political Science  Assoc. Director
1541 Lilac Lane, Room 504  Center for Research Methods
University of Kansas University of Kansas
http://pj.freefaculty.org   http://quant.ku.edu

[[alternative HTML version deleted]]

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

[R] [R-pkgs] update: rockchalk 1.8.0

[Paul Johnson]

This will appear on CRAN mirrors soon. It's my update for Spring, 2013. I keep
track of R problems that arise in the regression course and try to facilitate
them. There are functions for describing data, presenting regression plots and
tables, some regression diagnostics.

Most of the usages are illustrated in the vignette rockchalk.  That has plenty
of illustrations if you care to take a quick look.

If you did not try this before, here is why you might want to. The function that
the package was organized around originally was plotSlopes: draw predicted value
lines on the same scatterplot. Now this is generalized a bit, the moderator
variable can be numeric or factor and I have bent over backwards to make this
flexible for the end users. If you run the examples for predictOMatic and
plotSlopes and plotCurves, you will get the idea.

The rest is details.

I started a NEWS file using Emacs org mode, here it is:

* version 1.8

This is the end of the Spring semester, so its time for the new rockchalk
release.

** New, general, flexible framework for calculating marginal effects
   in regression models, linear or otherwise.

*** newdata function works. It can scan a regression, isolate the
   predictors, and then make a mix and match new data object for use
   with a predict function.  This is convenient for users but also very
   flexible.

*** The newdata framework is built on top of divider methods that can
check whether a variable is numeric or categorical, and select
example values according to user-specified criteria.

*** predictOMatic works dependably! Please try
example(predictOMatic). The problem with single predictor models
that bugged users of rockchalk 1.6.2 has been solved.

*** predictOMatic argument interval = c(none, confidence,
prediction).  Guess what that is supposed to do? For glm,
which does not provide a confidence interval, I've written code
for an approximate Wald type CI, and hope to do better in future.

** Regression diagnostics.

*** getPartialCor: get partial correlations from a fitted model
(student convenience).

*** getDeltaRsquare: Returns the change in estimated R-square observed
when each predictor is removed.  This is the squared semi-partial
correlation coefficients (student convenience).

*** mcDiangose for multicollinearity diagnostics (student convenience)

** MeanCenter: add arguments to make selection of variables for
centering more convenient when users don't like the automatic
options centerOnlyInteractors.

** plotSlopes, plotCurves:
 *** Added intervals argument, for confidence and prediction intervals.

*** Added opacity argument to determine darkness of interval regions
(which use the transparency alpha layer.).

*** A lot of fiddling under the hood to make colors consistent when
levels of modx are altered to compare plots of a given model.

*** Can produce a simple regression prediction plot if modx argument
is omitted. This is a widely requested feature.

Please run example(plotSlopes) and example(plotCurves)

*** Changes under the hood. The common plotting functions of
plotSlopes and plotCurves are abstracted into a function
plotFancy, so now this will be eaiser for me to maintain. The
plotting ritual is the same, why keep 2 functions, you ask?
plotCurves tolerates more complicated regression
formula. plotSlopes leads to testSlopes, and hence to
plot.testSlopes.

** addLines: communication between 2 dimensional regression plots and 3
dimensional plots from plotPlane. Run example(addLines).

** plot.testSlopes. Run testSlopes on an interactive model. For a
model with 2 ocontinuous predictors that interact, this will generate
an ABSOLUTELY EXCELLENT and highly informative plot displaying the
effect of the interaction.

** outreg: LaTeX tables from regression models.

*** Reworked with new arguments to make tables for more types
of regressions. There's quite a bit more room for users to customize
the type of diagnostics they want to report.

The wide variety of output types from regression models is very
bothersome. I refuse to write a separate outreg method for each
different regression packages. If you want to use a package
from an author who is willing to do that, consider the texreg
package.

*** outreg2HTML. converts outreg to HTML markup and into a file.
Intended for importation into word processor programs.

** New vignette Rstyle. Most of the source-code files have been reformatted
to comply with my advice.

** genCorrelatedData2

*** genCorrelatedData2. For regression examples, suppose you want to
have 6 columns of an MVN with a certain mean and covariance
structure. And you want the regression formula to have
interactions and squared-terms. No more hassle. This is a
framework that works. Users set the mean, standard deviations, and
correlation values in various ways. Run
example(genCorrelatedData2).

*** To support that, there are more generally useful

Re: [R] complexity of operations in R

[Paul Johnson]

On Thu, Jul 19, 2012 at 11:11 AM, Bert Gunter gunter.ber...@gene.com wrote:
 Hadley et. al:

 Indeed. And using a loop is a poor way to do it anyway.

 v - as.list(rep(FALSE,dotot))

 is way faster.

 -- Bert


Its not entirely clear to me what we are supposed to conclude about this.

I can confirm Bert's claim that his way is much  faster.  I ran
Johan's code, then Bert's suggestion.  The speedup is almost
unbelievable.


 system.time(f(dotot))

   user  system elapsed
 25.249   0.332  25.606
 #system.time(g(dotot))
 system.time(h(dotot))
   user  system elapsed
 15.753   0.404  16.173
 p - function(dotot){v - as.list(rep(FALSE,dotot))}
 system.time(p(dotot))
   user  system elapsed
  0.004   0.000   0.004


Why is this so much faster?


 f
function(dotot){
  v-matrix(ncol=1,nrow=0)
  for(i in 1:dotot){
v-rbind(v,FALSE)
  }
  return(v)
}

Obvious problems.

1. The return() at the end is unnecessary. Cutting that off saves 4%
or so on time.


 f2 - function(dotot){
+   v-matrix(ncol=1,nrow=0)
+   for(i in 1:dotot){
+ v-rbind(v,FALSE)
+   }
+   v
+ }
 system.time(f2(dotot))
   user  system elapsed
 24.150   0.220  24.391

2. Use of rbind is SLOW.  Calling rbind over and over again is VERY
slow (http://pj.freefaculty.org/R/WorkingExamples/stackListItems.Rout).

3. Calling matrix over and over to create nothing is SLOW.


If we just want a giant empty list, do this:

mylst - vector(mode=list, length=dotot)

 system.time(mylst - vector(mode=list, length=dotot))
   user  system elapsed
  0.000   0.000   0.001

In a practical application, where something has to be created to go
into the list, I think some speed considerations should be:

1. When a 'thing' is created and put in the list, do we need to
allocate memory twice like this

x - do.whatever()

mylst[[999]] - x

This will create x, and then make a duplicate of x to go into mylst. BORING!

if x is a big thing, I think better to avoid copying by making it anonymous:

mylst[[999]] - do.whatever()

If memory is a consideration at all, this is better, and it avoids the
problem of allocating memory twice.


Anyway, it seems to me the point of this thread should be that
allocating a big giant list of nothings is arbitrarily fast, but it is
still not known

1.  what is the efficient way to create things to go into the list,

2. what is the best way to make use of the results once they are collected.

I'm sure for looping and calling rbind lots of times is a really bad
idea.  Everything else is on the table for discussion. Its not the for
loop that makes the original f function slow, its repeated use of
matrix and rbind.

pj
-- 
Paul E. Johnson
Professor, Political ScienceAssoc. Director
1541 Lilac Lane, Room 504 Center for Research Methods
University of Kansas   University of Kansas
http://pj.freefaculty.orghttp://quant.ku.edu

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

[R] device mismatch, coordinates trouble with X11 and pdf devices

[Paul Johnson]

Greetings.

I'm trying to understand a problem on a Dell Laptop.  Details below,
also uploaded the R working example that I pasted below.

http://pj.freefaculty.org/scraps/testSymbols.R

 sessionInfo()
R version 2.15.1 (2012-06-22)
Platform: x86_64-pc-linux-gnu (64-bit)

locale:
 [1] LC_CTYPE=en_US.UTF-8   LC_NUMERIC=C
 [3] LC_TIME=en_US.UTF-8LC_COLLATE=en_US.UTF-8
 [5] LC_MONETARY=en_US.UTF-8LC_MESSAGES=en_US.UTF-8
 [7] LC_PAPER=C LC_NAME=C
 [9] LC_ADDRESS=C   LC_TELEPHONE=C
[11] LC_MEASUREMENT=en_US.UTF-8 LC_IDENTIFICATION=C

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


## Paul Johnson paulj...@ku.edu
## 2012-08-10
##
## I've got trouble with a mis-match between screen and pdf devices.
## Please run this and tell me if the point z's marker is on the
## intersection of the 2 circles for you.

## Why do I ask? On my Dell M4600 laptop, the marker (the pch=12) is
## on the intersection in the x11 device. See the screenshot:
##
## http://pj.freefaculty.org/scraps/testSymbol-01.png
##
## But the pdf output of same is not
##
## http://pj.freefaculty.org/scraps/testSymbol.pdf
##
## Notice one other weird thing. The circle sizes change. Maybe
## that's just the same weird thing.

## The X11 device is 7in x 7in. Right? Same as pdf.
## Why would coordinates from locator work in one display, but not the
## other?  I've abused the symbols function somehow? locator is device
## dependant?

## Maybe its the video driver  LCD.
## This is a Dell Precision M4600. The LCD display is 13.6 inches x
## 7.6 inches (345.4mm x 193mm). By my reckoning, it is not precisely
## 16:9.  The Nvidia card says the native resolution is 1920x1080,
## which is exactly 16:9.




saveFile - FALSE
if(saveFile){
pdf(file=testSymbol.pdf, onefile=FALSE, paper=special, height = 7,
width = 7)
}

plot(c(0,1), c(0,1), type = n, xlab = , ylab = )


x - c(0.3, 0.3)
y - c(0.7, 0.7)
points(x[1], x[2])
points(y[1], y[2])

symbols(x[1], x[2], circles=c(0.3), lty=2, inches = FALSE, add=TRUE)

symbols(y[1], y[2], circles=c(0.3), lty=2, inches = FALSE, add=TRUE)

##Here's what I get from locator(1) on the intersection
z - list(x=0.4028, y=0.6485)
if (is.null(z)) z - locator(1) ## click on intersection of lines


points(z, pch=12)
text(x =0.96*z[[1]], y = 1.05*z[[2]], label=z)

dev.off()




-- 
Paul E. Johnson
Professor, Political ScienceAssoc. Director
1541 Lilac Lane, Room 504 Center for Research Methods
University of Kansas   University of Kansas
http://pj.freefaculty.orghttp://quant.ku.edu

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

[R] shade overlapping portions of circles (or other shapes)

[Paul Johnson]

I'm making some illustrations and it would be convenient to
automatically shade the overlapping portions of circles.  These
illustrations are for Social Choice theory, a field in political
science and economics.  I've wrestled together some examples so you
can see what I mean, but have not mastered the color overlapping
sections problem (as you will see):

http://pj.freefaculty.org/scraps/t-m15g.pdf

The shaded area does not overlap perfectly, mainly because I isolated
the points manually with locator(). Tedious!

I've not found R functions just for this.  Are there any packages?

You can get an idea of the problem if you run, for example

library(plotrix)
example(draw.circle)

See the red and the blue circles overlap?  What's the best way to
shade in that overlapping area.

I don't have Mathematica anymore, but they seem to be able to do this
more-or-less easily

http://mathematica.stackexchange.com/questions/2554/how-to-plot-venn-diagrams-with-mathematica

I want to make R do same ( or develop package to do it).

I was under some pressure to get a graph finished that I linked to
above. But I would not want to do it that way again. The shapes won't
always  be perfect circles, they may be ellipses or other convex
shapes.  But if I could understand the problem for circles, maybe I
could solve for others.

Do you have suggestions?

I am considering this. Collect the co-ordiates from the circle points

 C1 - draw.circle(3.5,7,0.8,border=blue,lty=2,lwd=2)
 C2 - draw.circle(2.5,8,0.6,border=red,lty=3,lwd=3)

Can I isolate the relevant arcs from C1 and C2? If I could, I'd join them

Cjoin - list(x = c(C1$x, C2$x), y = c(C1$y, C2$y))

And use that for plotting.

pj


-- Paul E. Johnson
Professor, Political ScienceAssoc. Director
1541 Lilac Lane, Room 504 Center for Research Methods
University of Kansas   University of Kansas
http://pj.freefaculty.orghttp://quant.ku.edu

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

[R] premature evaluation of symbols. Is that the way to describe this problem?

[Paul Johnson]

Dear eveRybody

In the package rockchalk, I have functions that take regressions and
make tables and plots.  Sometimes I'll combine a model with various
arguments and pass the resulting list around to be executed, usually
using do.call.  While debugging a function on a particularly large
dataset, I noticed that I've caused inconvenience for myself. I get
the correct calculations, but the arguments into functions are not
symbols when the function call happens. They are fully written out.

So the debugger does not see function(x), it sees
function(c(1,3,1,3,45,2,4,2).  Debugging begins with a
comprehensive listing of the elements in every object, which freezes
Emacs/ESS and generally ruins my quality of life.

I don't know the right words for this, but it seems to me object names
are parsed and evaluated before they need to be. Function calls
include the evaluated arguments in them, not just the symbols for
them.  I have some guesses on fixes, wonder what you think. And I
wonder if fixing this problem might generally make my functions faster
and more efficient because I'm not passing gigantic collections of
numbers through the function call.  I know, I don't have all the right
words here.

I've built a toy example that will illustrate the problem, you tell me
the words for it.

## Paul Johnson 2014-04-10
dat - data.frame(x = rnorm(50),y = rnorm(50))
m1 - lm(y ~ x, dat)

myRegFit - function(model, nd) predict(model, nd)

mySpecialFeature - function(model, ci){
pargs - list(model, model.frame(model)[1:3, ])
res - do.call(myRegFit, pargs)
print(res)
}

mySpecialFeature (m1)

debug(myRegFit)

mySpecialFeature (m1)

Note when the debugger enters, it has the whole model's structure
splatted into the function call.

 dat - data.frame(x = rnorm(50),y = rnorm(50))
 m1 - lm(y ~ x, dat)
 myRegFit -function(model, nd) predict(model, nd)
 mySpecialFeature - function(model, ci){
+ pargs - list(model, model.frame(model)[1:3, ])
+ res - do.call(myRegFit, pargs)
+ print(res)
+ }
 mySpecialFeature (m1)
  1   2   3
-0.04755431  0.35162844 -0.11715522
 debug(myRegFit)
 mySpecialFeature (m1)
debugging in: (function (model, nd)
predict(model, nd))(list(coefficients = c(0.0636305741709566,
-0.177786836929453), residuals = c(-0.0152803151982162, -0.885875162659858,
-1.23645319405006, -1.77900639571943, -1.9952045397527, 1.38150266407176,
-2.27403449262599, 0.0367524776530579, -0.881037818467492, -1.10816713568432,
-0.55749829201921, -0.372526253742828, -0.353208893775679, 0.531708523456415,
-0.43187124865558, 1.03973431972897, 0.849170115617157, 1.11227803262189,
0.47216440383252, 0.920060697785203, -0.374672861268964, 2.94683565121636,
0.514112041811711, -0.52321362055969, -0.0412387814196237, 0.983863448669766,
0.534230127442599, -0.869960511196742, 1.90586406082412, -1.84705932449576,
0.806425475391075, 1.90939977897903, 0.41030042787483, 0.994503041407507,
0.715719209301158, -0.538096591457249, -0.482411304681239, 0.0323998214753804,
0.551162374882342, -0.618989357027834, 1.08996565055366, -0.697423620816604,
1.38170655971013, 1.55752893685726, -0.0929258405664267, -1.00210610433922,
-1.51879925258188, -1.57050250989563, -1.06868502360026, 0.458860605094578
), effects = c(-0.661274203468574, -1.07255577360914, -1.36123824096605,
-1.71875465303308, -1.8806486154155, 1.38636416232103, -2.29259163100096,
0.153263315278269, -0.950879523052079, -0.963705724647863, -0.62175976245114,
-0.423965256680951, -0.350662885659068, 0.469175412149025, -0.400505083448679,
1.03440116252973, 0.878739280288923, 1.16574672001397, 0.358222935858666,
1.00514946836967, -0.316592303881481, 2.8507611924072, 0.573209391002668,
-0.393720180068215, 0.0971873363200073, 1.23818281311352, 0.449576129222722,
-0.929511151618747, 1.97922180250824, -1.80820009744905, 0.877996855335966,
1.86871623414376, 0.226023354471842, 0.814892815951223, 0.821400980265047,
-0.536299037896556, -0.358703204255386, 0.105714598012197, 0.543301738010905,
-0.643659132172249, 1.26412624281219, -0.808498804261978, 1.32273956796476,
1.57655585529458, 0.022266343917185, -1.20958321975888, -1.52288584310647,
-1.60904879400386, -1.08384090772898, 0.567729611277337), rank = 2L,
fitted.values = c(-0.0475543078742839, 0.351628437239565,
-0.11715522155, 0.165041007767092, 0.247859323684996,
0.0805663573239046, 0.0448510221995737, 0.250840726246576,
-0.0333621333428176, 0.293467635417987, -0.0248518201238109,
-0.00529650906918994, 0.0770350456001708, -0.0222159311803766,
0.120988140963125, 0.0650186744394092, 0.118247568257686,
0.154696293986828, -0.100617877288265, 0.202919505525394,
0.161729772710581, -0.0733692278375946, 0.163280463324874,
0.270640256833884, 0.284263319806951, 0.461009994308065,
-0.0559520919143943, -0.0176674192887905, 0.185028727541523,
0.132415672762266, 0.182304379869478, 0.011106443703975,
-0.207885424899216, -0.200768100305762, 0.234325514404888

Re: [R] very basic HLM question

[Paul Johnson]

2011/2/5 Sebastián Daza sebastian.d...@gmail.com:
 Hi everyone,

 I need to get a between-component variance (e.g. random effects Anova), but
 using lmer I don't get the same results (variance component) than using
 random effects Anova. I am using a database of students, clustered on
 schools (there is not the same number of students by school).

 According to the ICC1 command, the interclass correlation is .44

 ICC1(anova1)
 [1] 0.4414491

If you don't tell us exactly what model you are calculating in
anova1, how would we guess if there is something wrong?

Similarly, I get this
 ICC1
Error: object 'ICC1' not found

so it must mean you've loaded a package or written a function, which
you've not shown us.

I googled my way to a package called multilevel that has ICC1, and
its code for ICC1 shows a formula that does not match the one you used
to calculate ICC from lmer.

function (object)
{
MOD - summary(object)
MSB - MOD[[1]][1, 3]
MSW - MOD[[1]][2, 3]
GSIZE - (MOD[[1]][2, 1] + (MOD[[1]][1, 1] + 1))/(MOD[[1]][1,
1] + 1)
OUT - (MSB - MSW)/(MSB + ((GSIZE - 1) * MSW))
return(OUT)
}

I'm not saying that's right or wrong, just not obviously identical to
the formula you proposed.


 However, I cannot get the same ICC from the lmer output:

 anova2 - lmer(math ~ 1 + (1|schoolid), data=nels88)
 summary(anova2 - lmer(math ~ 1 + (1|schoolid), data=nels88))


Instead, do this (same thing, fits model only once):

 anova2 - lmer(math ~ 1 + (1|schoolid), data=nels88)
 summary(anova2)

Note that lmer is going to estimate a normally distributed random
effect for each school, as well as an individual observation random
effect (usual error term) that is assumed independent of the
school-level effect.  What is anova1 estimating?


-- 
Paul E. Johnson
Professor, Political Science
1541 Lilac Lane, Room 504
University of Kansas

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

Re: [R] different results in MASS's mca and SAS's corresp

[Paul Johnson]

On Sat, Feb 5, 2011 at 9:19 AM, David Winsemius dwinsem...@comcast.net wrote:

 On Feb 4, 2011, at 7:06 PM, Gong-Yi Liao wrote:

 Dear list:

  I have tried MASS's mca function and SAS's PROC corresp on the
  farms data (included in MASS, also used as mca's example), the
  results are different:

  R: mca(farms)$rs:
             1             2
 1   0.059296637  0.0455871427
 2   0.043077902 -0.0354728795
 3   0.059834286  0.0730485572
 4   0.059834286  0.0730485572
[snip]

     And in SAS's corresp output:

                               Row Coordinates

                                       Dim1       Dim2

                         1           1.0607    -0.8155
                         2           0.7706     0.6346
                         3           1.0703    -1.3067
                         4           1.0703    -1.3067
                         5           0.2308     0.9000
[snip]
       Is MASS's mca developed with different definition to SAS's
       corresp ?

 No, it's just that the values can only be defined up to a scaling factor
 (the same situation as with eigenvector decompostion). Take a look at the
 two dimensions, when each is put on the same scale:

 cbind(scale(rmca$D1),scale(smca$Dim1) )
            [,1]        [,2]
  [1,]  1.2824421  1.28242560
  [2,]  0.9316703  0.93168561
  [3,]  1.2940701  1.29403231
  [4,]  1.2940701  1.29403231
  [5,]  0.2789996  0.27905048

 cbind(scale(rmca$D2),scale(smca$Dim2) )
             [,1]        [,2]
  [1,]  1.06673426 -1.06677626
  [2,] -0.83006158  0.83012474
  [3,]  1.70932841 -1.70932351
  [4,]  1.70932841 -1.70932351
  [5,] -1.17729983  1.17729909

 David Winsemius, MD
 West Hartford, CT


When I came to David's comment, I understood the theory, but not the
numbers in his answer.  I wanted to see the MASS mca answers match
up with SAS, and the example did not (yet).

Below see that if you scale the mca output, and then multiply column 1
of the scaled results by 0.827094, then  you DO reproduce the SAS
column 1 results exactly.  Just rescale item 1 in mca's first column
to match the SAS output.  Repeat same with column 2, multiply
-0.7644828, and you reproduce column 2.

 rmca - mca(farms)
 scalermca - scale(rmca$rs)
 scalermca[1,]
   12
1.282442 1.066734
 1.0607/1.282442
[1] 0.827094
 -0.8155/1.06673426
[1] -0.7644828
 cbind(scalermca[,1] * 0.827094, scalermca[,2] *  -0.7644828)
  [,1][,2]
1   1.06070017 -0.8154
2   0.77057891  0.63456780
3   1.07031764 -1.30675217
4   1.07031764 -1.30675217
5   0.23075886  0.90002547
6   0.6943  0.60993995
7   0.10530240  0.78445402
8  -0.27026650  0.44225049
9   0.13426089  1.15670532
10  0.11861965  0.64778456
11  0.23807570  1.21775202
12  1.01156703 -0.01927226
13  0.28051938 -0.59805897
14 -1.17343686 -0.27122981
15 -0.83838041 -0.64003061
16 -0.05453708 -0.22925816
17 -0.91732401 -0.49899374
18 -0.92694148 -0.00774156
19 -1.30251038 -0.34994509
20 -1.30251038 -0.34994509

So, that does reproduce SAS exactly.  And I'm a little frustrated I
can't remember the matrix command to get that multiplication done
without cbinding the 2 columns together that way.

Question: I don't use mca, but to people who do, how are results
supposed to be scaled?  Is there a community accepted method or is
every user on his/her own to fiddle up the numbers however?

-- 
Paul E. Johnson
Professor, Political Science
1541 Lilac Lane, Room 504
University of Kansas

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

Re: [R] boostrap an nls regression

[Paul Johnson]

Hello!

On Thu, Feb 3, 2011 at 6:27 AM, m234 mhairialexan...@gmail.com wrote:

 II functional response for 2 data sets:

 nls(eaten~(a*suppl)/(1+a*h*suppl)

 where eaten is the number of prey eaten by a predator and suppl is the
 number of prey initially supplied to the same predator.

 I have parameter estimates of 'a' and 'h' for the two populations studied
 and would like to know if there is a significant different in the estimates
 between the two i.e. a1 vs a2, h1 vs h2. I would like to bootstap the data
 to get multiple (~1000) estimates and compare then via a ttest or
 equivalent. Is it possible to do this and obtain multiple estimations which

Hi,

Please read the posting guide--a complete example of your code a link
to a data frame will get you more helpful answers.

I can tell you how to do bootstrap re-sampling on a model that is
estimated from one data set. I wrote up some notes on that last summer
(go to middle of this
http://pj.freefaculty.org/R/SummerCamp2010/PJ-Lectures/functions1.pdf).
I have other intro R material floating about under the R part of that
link.

But I'm a bit stumped by your question for statistical reasons--
nothing to do with R.  From a statistical point of view, how do you
test the difference between coefficients from 2 different fitted nls
models, which are based on separate data samples?  I don't think
that's a trivial stat question, even if you have large samples and you
only need to calculate that estimated difference one time.

Then I wonder, how would a person do re-sampling when there are 2 data
sets involved. If you know of a cite on that, I'd like to see it.

One avenue I'd consider is to stack the 2 data sets together and
rewrite my nls to estimate one equation with indicator variables
(dummy variables) to separate the estimates for the 2 separate
equations. But there would be some pooling tests you'd have to run
first, to justify the idea that the 2 sets of data belong in the same
analysis in the first place.

Know what I mean? Suppose eaten is one long column, for both sets
combined, but create suppl1 and suppl2 that are 0 for the other data
set's cases, but suppl for the right one.

Fit this:

combmod -  nls(eaten~(a1*suppl1)/(1+a1*h1*suppl1) +
(a2*suppl2)/(1+a2*h2*suppl2))

This would conceivably allow comparison of a1 and a2. I think. I'm
trying to remember sampling theory on nls.

Well, in summary, I think you've got a harder stat problem than you
have R problem.  If you write out the code you use for a whole
exercise to do this 1 time, we might see what to do.  But remember to
post the full working example--as much as you have, anyway.

-- 
Paul E. Johnson
Professor, Political Science
1541 Lilac Lane, Room 504
University of Kansas

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

[R] example to demonstrate benefits of poly in regression?

[Paul Johnson]

Here's my little discussion example for a quadratic regression:

http://pj.freefaculty.org/R/WorkingExamples/regression-quadratic-1.R

Students press me to know the benefits of poly() over the more obvious
regression formulas.

I think I understand the theory on why poly() should be more numerically
stable, but I'm having trouble writing down an example that proves the
benefit of this.

I am beginning to suspect that because R's lm uses a QR decomposition under
the hood, the instability that we used to see when we inverted (X'X) is no
longer so serious as it used to be. When the columns in X are

x1 x2 x3 x3squared

then the QR decomposition is doing about the same thing as we would do if
we replaced x3 and x3squared by the poly(x3, 2) results.

Or not. I'm not arguing that, I'm thinking out loud, hoping you'll correct
me.

pj


-- 
Paul E. Johnson
Professor, Political Science  Assoc. Director
1541 Lilac Lane, Room 504  Center for Research Methods
University of Kansas University of Kansas
http://pj.freefaculty.org   http://quant.ku.edu

[[alternative HTML version deleted]]

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

Re: [R] example to demonstrate benefits of poly in regression?

[Paul Johnson]

On Mon, Apr 1, 2013 at 12:42 PM, Gabor Grothendieck ggrothendi...@gmail.com
 wrote:

 On Mon, Apr 1, 2013 at 1:20 PM, Paul Johnson pauljoh...@gmail.com wrote:
  Here's my little discussion example for a quadratic regression:
 
  http://pj.freefaculty.org/R/WorkingExamples/regression-quadratic-1.R
 
  Students press me to know the benefits of poly() over the more obvious
  regression formulas.
 
  I think I understand the theory on why poly() should be more numerically
  stable, but I'm having trouble writing down an example that proves the
  benefit of this.
 
  I am beginning to suspect that because R's lm uses a QR decomposition
 under
  the hood, the instability that we used to see when we inverted (X'X) is
 no
  longer so serious as it used to be. When the columns in X are
 
  x1 x2 x3 x3squared
 
  then the QR decomposition is doing about the same thing as we would do if
  we replaced x3 and x3squared by the poly(x3, 2) results.
 
  Or not. I'm not arguing that, I'm thinking out loud, hoping you'll
 correct
  me.
 

 # 1. One benefit is if you fit a higher degree polynomial using
 poly(x, n) the lower degree coefficients will agree with the fit using
 a lower n.
 # For example, compare these two:
 set.seed(123)
 y - rnorm(10); x - 1:10
 lm(y ~ poly(x, 2))
 lm(y ~ poly(x, 3))

 # however, that is not true if you don't use orthogonal polynomials.
 Compare these two:
 lm(y ~ poly(x, 2, raw = TRUE))
 lm(y ~ poly(x, 3, raw = TRUE))


Thanks, Gabor:

As usual, you are very helpful. I can't thank you enough.

I'm pasting your comments to the end of that R example file I mentioned
before, changing the examples to fit with the variable names I use there.

I think your reason #2 is the prize winner, I can't see anything wrong with
that one.

Reason #1 is not quite as strong, only holds strictly if the poly() is the
only predictor.  Users usually have other variables in the model.

If you--or anybody else--has other ideas about this, I am delighted to hear
your ideas.  I still hold out hope that the numerical stability argument
will find some traction.

pj



 # 2. A second advantage is that the t statistics are invariant under
 shift if you use orthogonal polynomials
 # compare these two:
 summary(lm(y ~ poly(x, 2)))
 summary(lm(y ~ poly(x+1, 2)))

 # however, that is not true if you don;t use orthogonal polynomials,
 Compare these two:
 summary(lm(y ~ poly(x, 2, raw = TRUE)))
 summary(lm(y ~ poly(x+1, 2, raw = TRUE)))

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




-- 
Paul E. Johnson
Professor, Political Science  Assoc. Director
1541 Lilac Lane, Room 504  Center for Research Methods
University of Kansas University of Kansas
http://pj.freefaculty.org   http://quant.ku.edu

[[alternative HTML version deleted]]

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

Re: [R] DUD (Does not Use Derivatives) for nonlinear regression in R?

[Paul Johnson]

On Apr 1, 2013 1:10 AM, qi A send2...@gmail.com wrote:

 Hi, All

 SAS has DUD (Does not Use Derivatives)/Secant Method for nonlinear
 regression, does R offer this option for nonlinear regression?

 I have read the helpfile for nls() and could not find such option, any
 suggestion?


nelder-mead is default algorithm in optim. It does not use derivatives. dud
is from same generation, but John Nash recommended N-M method.

Pj

 Thanks,

 Derek

 [[alternative HTML version deleted]]

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

[[alternative HTML version deleted]]

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

Re: [R] Superscript

[Paul Johnson]

On Thu, Apr 4, 2013 at 4:39 AM, Shane Carey careys...@gmail.com wrote:

 Hi William,

 for (i in one:length(DATA_names))
   if ((grepl(_,DATA_names[i]))==TRUE)
 DATA_names[i]-f(DATA_names[i]))

 I keep getting an error saying: incompatible types (from symbol to
 character) in subassignment type fix


I would urge you not to obliterate DATA_names by writing on top of itself.
That creates bugs that are really tough to solve.

Remember that a character variable is not the same thing as an expression.
I mean, an expression is generally a list structure, not interchangeable
with a character string.  If that f() function is creating expressions, and
you are trying to obliterate a character vector with expressions, you will
get errors, as you have seen.

I'm pretty sure that to catch the results of a function that creates
expressions, your receiver object has to be a list.  My first try would be

myList - vector(list, length(DATA_names))

for (i in one:length(DATA_names))
  if ((grepl(_, DATA_names[i])) == TRUE)
myList[[i]] - f(DATA_names[i]))

I copied your code to my system and found I couldn't run it because you had
one more paren than was needed and the word one was undefined. But after
fixing that, here's evidence my way works:

Note I inserted some spaces in your code, since it is inhumane to smash
together all those letters around a poor little - :).


DATA_names-c(A_mgkg,B_mgkg,C_mgkg,D_mgkg,E_mgkg,F_mgkg,G_mgkg,H
mgkg)

f - function (name)
{
  # add other suffices and their corresponding plotmath expressions to the
list
  env - list2env(list(mgkg = bquote(mg ~ kg^{-1}),
   ugkg = bquote(mu * g ~ kg^{-1})),
  parent = emptyenv())
  pattern - paste0(_(, paste(objects(env), collapse=|), ))
  bquoteExpr - parse(text=gsub(pattern,
~(.(\\1)),
name))[[1]]
## I use do.call() to work around the fact that bquote's first argument is
## not evaluated.
  do.call(bquote, list(bquoteExpr, env))
}

myList - vector(list, length(DATA_names))

for (i in 1:length(DATA_names))
  if ((grepl(_, DATA_names[i])) == TRUE)
myList[[i]] - f(DATA_names[i])


myList


 myList
[[1]]
A ~ (mg ~ kg^{
-1
})

[[2]]
B ~ (mg ~ kg^{
-1
})

[[3]]
C ~ (mg ~ kg^{
-1
})

[[4]]
D ~ (mg ~ kg^{
-1
})

[[5]]
E ~ (mg ~ kg^{
-1
})

[[6]]
F ~ (mg ~ kg^{
-1
})

[[7]]
G ~ (mg ~ kg^{
-1
})

[[8]]
NULL


I asked a very similar question here last year and got many interesting
suggestions. I was so fascinated I wrote a little programming vignette
called Rchaeology in my package rockchalk. As soon as you settle on the
most direct route from start to end, please post, I might learn new tricks.

In your case, I've been wondering if you ought to edit the string before it
becomes an expression, or whether you should edit the expression after. I
needed a correct expression to start, though.  This writes the correct
thing in the middle of the graph. Correct?

xcorrect - expression(X ~~ mg kg^{-1})
plot(1:10, 1:10, type =n)
text(5, 5, xcorrect)

Supposing that's correct, watch this.

mgkg - quote(mg kg^{-1})

 text(4, 2, bquote(A ~~ .(mgkg)))
 text(4, 3, bquote(B ~~ .(mgkg)))
 text(4, 6, bquote(C ~~ .(mgkg)))

I've got the units captured in mgkg, and then when I want to draw that
into
and expression, i used bquote. substitute works as well, maybe more clear.

text(3, 8, substitute(B ~~ myUnit, list(myUnit = mgkg)))

In any case, I thought this was a fun question.

pj



 Have you any ideas on how to get around this, thanks again for your help,
 much appreciated.
 Cheers


 On Wed, Apr 3, 2013 at 5:33 PM, William Dunlap wdun...@tibco.com wrote:

  Are you trying to convert a column name like Na_mgkg to a plot label
  like Na (mg kg^-1) ?
  If so you will have to use both string manipulation functions like gsub()
  and expression manipulating
  functions like bquote().  E.g.,
 
  f - function (name)
  {
 # add other suffices and their corresponding plotmath expressions to
  the list
 env - list2env(list(mgkg = bquote(mg ~ kg^{-1}),
  ugkg = bquote(mu * g ~ kg^{-1})),
 parent = emptyenv())
 pattern - paste0(_(, paste(objects(env), collapse=|), ))
 bquoteExpr - parse(text=gsub(pattern,
   ~(.(\\1)),
   name))[[1]]
 # I use do.call() to work around the fact that bquote's first argument
  is not evaluated.
 do.call(bquote, list(bquoteExpr, env))
  }
 
  d - data.frame(Na_mgkg=1:10, K_ugkg=10:1)
  plot(Na_mgkg ~ K_ugkg, data=d, xlab=f(K_ugkg), ylab=f(Na_mgkg))
 
  Bill Dunlap
  Spotfire, TIBCO Software
  wdunlap tibco.com
 
 
   -Original Message-
   From: r-help-boun...@r-project.org [mailto:
 r-help-boun...@r-project.org]
  On Behalf
   Of Shane Carey
   Sent: Wednesday, April 03, 2013 8:02 AM
   To: r-help@r-project.org
   Subject: [R] Superscript
  
   Hi,
   How do I write a superscript within 

Re: [R] Make barplot with error bars, anova out of a table

[Paul Johnson]

Hi

On Thu, Apr 11, 2013 at 6:43 AM, Simza alex_steckba...@gmx.at wrote:
 Helo everybody,
 I'm new to R and have some issues with my data in R.

 My raw data look like that:

 ID Day size 1 1 7 1 1 7.2 1 1 7.1 2 1 7.3 2 1 7.4 2 1 7.2 3 1 7 3 1 7.1 3 1
 7.5 4 1 7.3 4 1 7.2 4 1 7.6 1 2 7 1 2 7.2 1 2 7.1 2 2 7.1 2 2 7.4 2 2 7.2 3
 2 7.5 3 2 7.1 3 2 7.5 4 2 7.2 4 2 7.2 4 2 7.3 1 3 7.4 1 3 7.2 1 3 7.1 2 3
 7.2 2 3 7.4 2 3 7.2 3 3 7.4 3 3 7.2 3 3 7.5 4 3 7.4 4 3 7.2 4 3 7.7

 What I want to do is:
 1) calculate the significant difference of the size for the 4 groups but
 each day separately!)? E.g. difference in size between ID 1:2, 1:3, 1:4,
 2:3, 2:4, 3:4 for each day.

 I already tried lm() and anova() but just get the difference between day, ID
 and size. I also tried to separate the data per hand but I think there
 should be a way in R to do that automatically!? Moreover, I was searching
 the web for similar examples but couln't find anything useful so far.

 2) to make a barplot with error bars of the standart error (from the mean).
 So far I used:
 barplot(matrix(c(Rtest.dat$pH.mean),nr=3), beside=T,
 col=c(black,grey,white), main=pH, names.arg=c(Green, Yellow,
 Blue, Red), ylab=pH) legend(topright, c(Day 1,Day 2,Day 3),
 cex=0.6, bty=n, fill=c(black,grey,white))
 But I have problems to add the error bars. Also here I already searched the
 web but couldn't manage to get a code working for my data.

R has a function called arrows that we use to draw intervals like
that. You have to calculate the intervals, of course, that's up to
you, and then draw them on the plot.

As luck would have it, I've been working a vaguely similar type of plot,

http://pj.freefaculty.org/scraps/plot-factor-proto.pdf

and I'm going to upload the R code that produces that plot. It shows
how to do arrows and such. I'm calling this a WorkingExample

http://pj.freefaculty.org/R/WorkingExamples/plot-regression-factors-1.R

 3) and the last thing I would need is to add * or letters (a,b,...) to the
 graph, indicating the significant difference between them.

Try text() to insert stuff on plots. * is writable by text.


 Hope someone can help me!

 Thanks for your help!!



 --
 View this message in context: 
 http://r.789695.n4.nabble.com/Make-barplot-with-error-bars-anova-out-of-a-table-tp4663979.html
 Sent from the R help mailing list archive at Nabble.com.

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



--
Paul E. Johnson
Professor, Political Science  Assoc. Director
1541 Lilac Lane, Room 504  Center for Research Methods
University of Kansas University of Kansas
http://pj.freefaculty.org   http://quant.ku.edu

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

[R] Needed: Beta Testers and Summer Camp Students

[Paul Johnson]

Greetings.

I'm teaching linear regression this semester and that means I write
more functions for my regression support package rockchalk.  I'm at
a point now were some fresh eyes would help, so if you are a student
in a course on regression, please consider looking over my package
overview document here:

http://pj.freefaculty.org/R/rockchalk.pdf

That tells how you can grab the test verion, which is now at 1.7.I'm
leaving rockchalk on CRAN at 1.6.3, but as the document linked above
explains, you can download the test version from our local repository
called kran. I have been making new packages every 3 days or so. If
you are a github user, you can clone a copy of the source if you like
(http://

The functions that have gotten the biggest workover are predictOMatic,
newdata, plotSlopes, plotCurves, and testSlopes. If you just install
the package and run those examples, you will be able to tell right
away if you are interested in trying to adapt these to your needs.
Generally speaking, I watch the students each semester to see which R
things frustrate them the most and then I try to automate them.
That's how the regression table function (outreg) started, and
difficulties in plotting predicted values for various kinds of
regression drive most of the rest. I guess the rockchalk vignette
explains all that.

If you are interested in that flavor of R, or know other people who
might be, spread the word:

we are offering a one-week summer course at the University of Kansas.

There is a discount for enrollment before the end of the month.  The R
course is part of our larger Summer Statistical Institute, which has
been growing rapidly for the past decade. We've had very popular
courses on structural equations modeling and hierarchical models.

Here's the more formal announcement.

Stats Camp 2013 registration is now in full swing. Last year we had
over 300 participants attend 11 different courses. This coming June we
have 15 different courses being offered. Please visit our web pages at
http://crmda.ku.edu/statscamp for information on the courses, a brief
syllabus for each, and registration information.

pj

--
Paul E. Johnson
Professor, Political Science  Assoc. Director
1541 Lilac Lane, Room 504  Center for Research Methods
University of Kansas University of Kansas
http://pj.freefaculty.org   http://quant.ku.edu

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

Re: [R] pglm package: fitted values and residuals

[Paul Johnson]

On Wed, Apr 24, 2013 at 3:11 AM,  alfonso.carf...@uniparthenope.it wrote:
 I'm using the package pglm and I'have estimated a random probit model.
 I need to save in a vector the fitted values and the residuals of the model
 but I can not do it.

 I tried with the command fitted.values using the following procedure without
 results:

This is one of those ask the pglm authors questions. You should take
it up with the authors of the package.  There is a specialized email
list R-sig-mixed where you will find more people working on this exact
same thing.

pglm looks like fun to me, but it is not quite done, so far as I can
tell. Well, the authors have not gone the extra step to make their
regression objects behave like other R regression objects.   In case
you need alternative software, ask in R-sig-mixed. You'll learn that
most of these can be estimated with other packages. But I really like
the homogeneous user interface that is spelled out in pglm, and I
expect my students will run into the same questions that you have..

I just downloaded their source code, you probably ought to do that so
you can understand what they are doing.   They provide the fitting
functions, but they do not do any of the other work necessary to make
these functions fit together with the R class framework.  There are no
methods for predict, anova, and so forth.

I'm in their R folder looking for implementations:

pauljohn@pols110:/tmp/pglm/R$ grep summary *
pauljohn@pols110:/tmp/pglm/R$ grep predict *
pauljohn@pols110:/tmp/pglm/R$ grep class *
pauljohn@pols110:/tmp/pglm/R$ grep fitted *
pglm.R:  # glm models can be fitted

Run

 example(pglm)

what can we do after that?

 plot(anb)
Error in xy.coords(x, y, xlabel, ylabel, log) :
  'x' is a list, but does not have components 'x' and 'y'

## Nothing.
## We do get a regression summary object, that's better than some
packages provide:

 anbsum - summary(anb)

## And a coefficient table

 coef(anbsum)
 Estimate  Std. error   t value  Pr( t)
(Intercept) -6.933764e-01 0.061391429 -11.294351205 1.399336e-29
wage 1.517009e-02 0.006375966   2.379261231 1.734738e-02
exper1.314229e-03 0.007400129   0.177595444 8.590407e-01
ruralyes-8.594328e-05 0.051334716  -0.001674175 9.986642e-01

 model.matrix(anb)
Error in terms.default(object) : no terms component nor attribute
 anova(anb)
Error in UseMethod(anova) :
  no applicable method for 'anova' applied to an object of class
c('maxLik', 'maxim')
 predict(anb)
Error in UseMethod(predict) :
  no applicable method for 'predict' applied to an object of class
c('maxLik', 'maxim')

So, if you want those features with these models, you'll have to get
busy and do a lot of coding!

While working on regression support lately, I've reached the
conclusion that if an R package that claims to do regression but
does not provide methods for summary, predict, anova, nobs, fitted,
logLik, AIC, and so forth, then it is not done yet. Otherwise, users
like you who expect to be able to run methods like fitted or such have
a bad experience, as you are having now.

Maybe somebody reading this will remind us where the common list of R
regression methods is listed. I know for sure I've seen a document
about these things, but I'm baffled now trying to find it. But I'm
sure there is one.


pj

 library(pglm)

 m1_S-pglm(Feed ~ Cons_PC_1 + imp_gen_1 + LGDP_PC_1 + lnEI_1 +
 SH_Ren_1,data,family=binomial(probit),model=random,method=bfgs,index=c(Year,IDCountry))

 m1_S$fitted.values
 residuals(m1)


 Can someone help me about it?

 Thanks


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

Re: [R] pglm package: fitted values and residuals

[Paul Johnson]

On Wed, Apr 24, 2013 at 4:37 PM, Achim Zeileis achim.zeil...@uibk.ac.at wrote:
 On Wed, 24 Apr 2013, Paul Johnson wrote:

 On Wed, Apr 24, 2013 at 3:11 AM,  alfonso.carf...@uniparthenope.it
 wrote:

 I'm using the package pglm and I'have estimated a random probit model.
 I need to save in a vector the fitted values and the residuals of the model
 but I can not do it.

 I tried with the command fitted.values using the following procedure
 without results:

 This is one of those ask the pglm authors questions. You should take it
 up with the authors of the package.  There is a specialized email list
 R-sig-mixed where you will find more people working on this exact same
 thing.

 pglm looks like fun to me, but it is not quite done, so far as I can tell.

 I'm sure that there are many. One of my attempts to write up a list is in
 Table 1 of vignette(betareg, package = betareg).

Yes! That's exactly the list I was thinking of.  It was driving me
crazy I could not find it.

Thanks for the explanation.  I don't think I should have implied that
the pglm author must actually implement all the methods, it is
certainly acceptable to leverage the methods that exist.  It just
happened that the ones I tested were not implemented by any of the
affiliated packages.

But this thread leads me to one question I've wondered about recently.

Suppose I run somebody's regression function and out comes an object.

Do we have a way to ask that object what are all of the methods that
might apply to you?  Here's why I wondered. You've noticed that
predict.lm has the interval=confidence argument, but predict.glm
does not. So if I receive a regression model, I'd like to say to it
do you have a predict method and if I could get that predict method,
I could check to see if there is a formal argument interval. If it
does not, maybe I'd craft one for them.

pj



 Personally, I don't write anova() methods for my model objects because I can
 leverage lrtest() and waldtest() from lmtest and linearHypothesis() and
 deltaMethod() from car as long as certain standard methods are available,
 including coef(), vcov(), logLik(), etc.

 Similarly, an AIC() method is typically not needed as long as logLik() is
 available. And BIC() works if nobs() is available in addition.

 Best,
 Z


 pj

 library(pglm)

 m1_S-pglm(Feed ~ Cons_PC_1 + imp_gen_1 + LGDP_PC_1 + lnEI_1 +

 SH_Ren_1,data,family=binomial(probit),model=random,method=bfgs,index=c(Year,IDCountry))

 m1_S$fitted.values
 residuals(m1)


 Can someone help me about it?

 Thanks


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





--
Paul E. Johnson
Professor, Political Science  Assoc. Director
1541 Lilac Lane, Room 504  Center for Research Methods
University of Kansas University of Kansas
http://pj.freefaculty.org   http://quant.ku.edu

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

[R] Trouble with methods() after loading gdata package.

[Paul Johnson]

Greetings to r-help land.

I've run into some program crashes and I've traced them back to methods()
behavior
after the package gdata is loaded.  I provide now a minimal re-producible
example. This seems bugish to me. How about you?

dat - data.frame(x = rnorm(100), y = rnorm(100))
lm1 - lm(y ~ x, data = dat)

methods(class = lm)

## OK so far

library(gdata)
methods(class = lm)
## epic fail



## OUTPUT.

 dat - data.frame(x = rnorm(100), y = rnorm(100))
 lm1 - lm(y ~ x, data = dat)

 methods(class = lm)
 [1] add1.lm*   alias.lm*  anova.lm   case.names.lm*
 [5] confint.lm*cooks.distance.lm* deviance.lm*   dfbeta.lm*
 [9] dfbetas.lm*drop1.lm*  dummy.coef.lm* effects.lm*
[13] extractAIC.lm* family.lm* formula.lm*hatvalues.lm
[17] influence.lm*  kappa.lm   labels.lm* logLik.lm*
[21] model.frame.lm model.matrix.lmnobs.lm*   plot.lm
[25] predict.lm print.lm   proj.lm*   qr.lm*
[29] residuals.lm   rstandard.lm   rstudent.lmsimulate.lm*
[33] summary.lm variable.names.lm* vcov.lm*

   Non-visible functions are asterisked

 library(gdata)
gdata: read.xls support for 'XLS' (Excel 97-2004) files ENABLED.

gdata: read.xls support for 'XLSX' (Excel 2007+) files ENABLED.

Attaching package: ‘gdata’

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

nobs

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

object.size

 methods(class = lm)
Error in data.frame(visible = rep.int(FALSE, n2), from = rep.int(msg,  :
  duplicate row.names: nobs.lm

 sessionInfo()
R version 3.0.0 (2013-04-03)
Platform: x86_64-pc-linux-gnu (64-bit)

locale:
 [1] LC_CTYPE=en_US.UTF-8   LC_NUMERIC=C
 [3] LC_TIME=en_US.UTF-8LC_COLLATE=en_US.UTF-8
 [5] LC_MONETARY=en_US.UTF-8LC_MESSAGES=en_US.UTF-8
 [7] LC_PAPER=C LC_NAME=C
 [9] LC_ADDRESS=C   LC_TELEPHONE=C
[11] LC_MEASUREMENT=en_US.UTF-8 LC_IDENTIFICATION=C

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

other attached packages:
[1] gdata_2.12.0.2

loaded via a namespace (and not attached):
[1] gtools_2.7.1 tcltk_3.0.0  tools_3.0.0


gdata is one of my favorite packages, its worth the effort to get to the
bottom of this.

--
Paul E. Johnson
Professor, Political Science  Assoc. Director
1541 Lilac Lane, Room 504  Center for Research Methods
University of Kansas University of Kansas
http://pj.freefaculty.org   http://quant.ku.edu

[[alternative HTML version deleted]]

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

Re: [R] significantly different from one (not zero) using lm

[Paul Johnson]

It is easy to construct your own test. I test against null of 0 first so I
can be sure I match the right result from summary.lm.

## get the standard error
seofb - sqrt(diag(vcov(lm1)))
## calculate t. Replace 0 by your null
myt - (coef(lm1) - 0)/seofb
mypval - 2*pt(abs(myt), lower.tail = FALSE, df = lm1$df.residual)

## Note you can pass a vector of different nulls for the coefficients
myt - (coef(lm1)  - c(0,1))/seofb

We could write this into a function if we wanted to get busy.  Not a bad
little homework exercise, I think.




 dat - data.frame(x = rnorm(100), y = rnorm(100))
 lm1 - lm(y ~ x, data = dat)
 summary(lm1)

Call:
lm(formula = y ~ x, data = dat)

Residuals:
Min  1Q  Median  3Q Max
-3.0696 -0.5833  0.1351  0.7162  2.3229

Coefficients:
 Estimate Std. Error t value Pr(|t|)
(Intercept) -0.001499   0.104865  -0.0140.989
x   -0.039324   0.113486  -0.3470.730

Residual standard error: 1.024 on 98 degrees of freedom
Multiple R-squared: 0.001224,Adjusted R-squared: -0.008968
F-statistic: 0.1201 on 1 and 98 DF,  p-value: 0.7297

 seofb - sqrt(diag(vcov(lm1)))
 myt - (coef(lm1) - 0)/seofb
 mypval - 2*pt(abs(myt), lower.tail = FALSE, df = lm1$df.residual)
 myt
(Intercept)   x
-0.01429604 -0.34650900
 mypval
(Intercept)   x
  0.9886229   0.7297031
 myt - (coef(lm1) - 1)/seofb
 mypval - 2*pt(abs(myt), lower.tail = FALSE, df = lm1$df.residual)
 myt
(Intercept)   x
  -9.550359   -9.158166
 mypval
 (Intercept)x
1.145542e-15 8.126553e-15


On Tue, Apr 30, 2013 at 9:07 PM, Elaine Kuo elaine.kuo...@gmail.com wrote:

 Hello,



 I am work with a linear regression model:

 y=ax+b with the function of lm.

 y= observed migration distance of butterflies

 x= predicted migration distance of butterflies



 Usually the result will show

 if the linear term a is significantly different from zero based on the
 p-value.

 Now I would like to test if the linear term is significantly different from
 one.

 (because I want to know if the regression line (y=ax+b) is significantly
 from the line with the linear term =1 and the intercept =0)



 Please kindly advise if it is possible

 to adjust some default parameters in the function to achieve the goal.

 Thank you.


 Elaine

 [[alternative HTML version deleted]]

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




-- 
Paul E. Johnson
Professor, Political Science  Assoc. Director
1541 Lilac Lane, Room 504  Center for Research Methods
University of Kansas University of Kansas
http://pj.freefaculty.org   http://quant.ku.edu

[[alternative HTML version deleted]]

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

[R] Automatic way to embed fonts into all pdf output?

[Paul Johnson]

This is a feature request. Or else a howto request.  Can there be some
simple, automatic way to make fonts embed into pdf output? Could this be in
options() or par()?

Why?

I recently wanted to create a paper for a conference in a format called
AAAI (which I had never heard of before because I live in a cave).

http://www.aaai.org/Publications/Author/author.php

In their instructions for LaTeX, it has the statement:

All fonts must be embedded in the PDF file — this includes your figures.

I remember that after one creates a pdf file in R with the pdf device, one
can then add fonts to it with the embedFonts() function. But I forget!  And
this is difficult to explain to nonprogramemrs, who just want to put some
figures in a document.   I wish there were something I could put in
options() or par() or somewhere to just make this always happen.  I
understand this will make pdf files larger, but I'm willing to spend the
disk.

-- 
Paul E. Johnson
Professor, Political Science  Assoc. Director
1541 Lilac Lane, Room 504  Center for Research Methods
University of Kansas University of Kansas
http://pj.freefaculty.org   http://quant.ku.edu

[[alternative HTML version deleted]]

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

Re: [R] R installation error

[Paul Johnson]

Meenu:

You have an elementary Linux setup and configuration problem to understand
first, before you worry about configuring and compiling your own R.  I
agree strongly that this is something that all linux users should learn to
do, but compiling R itself is like climbing Mt Everest as your first
mountain climb.

So, for whatever Linux us use, join their email support lists. Find out
more about how to install packages there. Not all Linux systems use yum
or apt, we can't help you with that.  SO learn about whatever linux you
have and find out if it separately installs packages required for compiling
programs.  It may even be you have a LInux system without sudo installed,
so learn how to install that.

After you succeed in compiling something simple, then get serious.

My course notes on this:

http://pj.freefaculty.org/guides/Computing-HOWTO/IntroTerminal-3/terminal-3.pdf




On Fri, May 10, 2013 at 7:33 AM, Meenu Chopra meenu.bioinf...@gmail.comwrote:

 Thanks to all

 The main problem with ma linux system is that its not able to install any
 software using sudo command.
 like i used command yum search libX11 , it shown that yum is not installed
 and when i use sudo apt-get install yum
 its giving error
 E: Unable to locate package yum
 same problem with all software.

 even when u use sudo apt-get update it shows error like

 W: Failed to fetch

 http://old.releases.ubuntu.com/ubuntu/dists/maverick-security/universe/binary-amd64/Packages.gz
 404  Not Found

 I am not getting whats wrong with my system

 Please help me out



 On Fri, May 10, 2013 at 2:32 AM, Jim Lemon j...@bitwrit.com.au wrote:

  On 05/09/2013 11:06 PM, Meenu Chopra wrote:
 
  Hiii
 
  I am trying to install R-2.15.2
  after doing ./configure its showing error: --with-x=yes (default) and
 X11
  headers/libs are not available
 
  and when i am running make its showing
  make: *** No targets specified and no makefile found. Stop.
  Even I read the install file also but i am not getting any solution.
 
   Hi Meenu,
  The second error is due to the first, so don't worry about it for the
  moment. Since you are compiling, you are probably on a Linux system. Try
  this:
 
  su
  root password
  yum search libX11
 
  this will produce a list of about 6 packages
  if you are on a 32 bit system
 
  yum install libX11-devel.i686 libX11.i686 libX11-common.noarch
 
  if you are on a 64 bit system
 
  yum install libX11-devel.x86_64 libX11.x86_64 libX11-common.noarch
 
  exit
 
  This should get the necessary stuff. You may also have problems with
  readline. If so, let us know.
 
  Jim
 



 --
 Meenu Chopra
 Research Associate
 Animal Genomics Lab
 NDRI, Karnal

 [[alternative HTML version deleted]]

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




-- 
Paul E. Johnson
Professor, Political Science  Assoc. Director
1541 Lilac Lane, Room 504  Center for Research Methods
University of Kansas University of Kansas
http://pj.freefaculty.org   http://quant.ku.edu

[[alternative HTML version deleted]]

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

[R] [R-pkgs] rockchalk_1.5.4 posted

[Paul Johnson]

Greetings:

rockchalk is a collection of functions to facilitate presentation of regression 
models.
It includes some functions that I have been circulating for quite some time 
(such as
outreg) as well as several others. The main aim is to allow people who do not
understand very much R to survive a course in intermediate regression analysis.
The examples included with the functions include more than the usual amount of 
detail.

This version features

1) a full vignette called rockchalk that illustrates many of the
functions in the package.  The vignette includes an explanation of why 
mean-centering
and residual-centering do not help with the inessential multicollinearity 
problem
that concerns users of regression models in which there are interactions or
squared terms.

2) a function summarize that is intended to remedy some of the shortcomings I
perceive in R's summary function.

The package is not coordinated with any particular textbook and should help in
any class in which students are expected to estimate regression models, 
summarize
them in professional-looking tables, and conduct various follow up diagnostics.

Highlights:

plotSlopes  testSlopes: draw simple slopes graphs and conduct hypothesis 
tests
for moderator variables.

plotCurves: an attempt to provide a termplot replacement that allows plots 
that
incorporate mediators.

plotPlane: a wrapper for persp to plot regressions in 3d (similar to 
scatterplot3d)

standardize, meanCenter, residualCenter: after fitting a non-centered model, 
use these
to experiment with the benefits (or lack thereof) that flow from centering.

-- 
Paul E. Johnson email: paulj...@ku.edu
http://pj.freefaculty.org   Assoc. Director
Professor, Political ScienceCenter for Research Methods and Data Analysis
1541 Lilac Lane, Rm 504 1425 Jayhawk Blvd.  
University of KansasWatson Library, Rm. 470 
Lawrence, Kansas 66045-3129 Lawrence, Kansas 66045-7555
Ph: (785) 864-3523  Ph: (785) 864-3353

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

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

[R] predictOMatic for regression. Please try and advise me

[Paul Johnson]

I'm pasting below a working R file featuring a function I'd like to polish up.

I'm teaching regression this semester and every time I come to
something that is very difficult to explain in class, I try to
simplify it by writing an R function (eventually into my package
rockchalk).  Students have a difficult time with predict and newdata
objects, so right now I'm concentrated on that problem.

Fitting regressions is pretty easy, but getting predicted values for
particular combinations of input values can be difficult when the
model is complicated. If you agree, please try out the following and
let me know how its use could be enhanced, or what other features you
might want.

I know folks are busy, so to save you the trouble of actually running
the code, I also paste in a session demonstrating one run through.

Here's the code:

##Paul Johnson
## 2012-04-20
## Facilitate creation of newdata objects for use in predict
## statements (for interpretation of regression models)

## newdataMaker creates the newdata frame required in predict.
## It supplies a data frame with all input variables at a
## central value (mean, mode) except for some variables that
## the user wants to focus on. User should
## supply a fitted model model and a focus list fl of variable
## values. See examples below.  The list must be a named list,
## using names of variables from the regression formula.
## It is not needed to call this
## directly if one is satisfied with the results from
## predictOMatic

newdataMaker - function (model = NULL, fl = NULL){
mf - model.frame(model)
mterms - terms(model)
mfnames - colnames(mf)
modelcv - centralValues(mf)
if (sum(!names(fl) %in% mfnames)  0) stop(cat(c(Error. The focus
list (fl) requests variables that are not included in the original
model. The names of the variables in the focus list be drawn from this
list: ,  mfnames)))
## Consider padding range of fl for numeric variables so that we
## get newdata objects including the min and max values.
mixAndMatch - expand.grid(fl)
mixAndMatch$combo - 1:nrow(mixAndMatch)
newdf - cbind(mixAndMatch, modelcv[  ,!colnames(modelcv) %in%
colnames(mixAndMatch)])
newdf
}



predictOMatic - function(model = NULL, fl = NULL, ...){
nd - newdataMaker(model, fl)
fit - predict(model, newdata=nd, ...)
cbind(fit, nd)
}



set.seed(12345)
x1 - rnorm(100)
x2 - rnorm(100)
x3 - rnorm(100)
x4 - rnorm(100)
x5 - rpois(100, 5)
x6 - rgamma(100, 2,1)
xcat1 - gl(2,50, labels=c(M,F))
xcat2 - cut(rnorm(100), breaks=c(-Inf, 0, 0.4, 0.9, 1, Inf),
labels=c(R, M, D, P, G))
dat - data.frame(x1, x2, x3, x4, x5, xcat1, xcat2)
rm(x1, x2, x3, x4, xcat1, xcat2)
dat$y - with(dat, 0.03 + 0.1*x1 + 0.1*x2 + 0.4*x3 -0.1*x4 + 0.01*x5 +
0.1*x6 + 2*rnorm(100))

##ordinary regression.
m1 - lm(y ~ x1 + x2 + x3 + x4 + x5 + x6 + xcat1 + xcat2, data=dat)
summary(m1)
## Use the fitted model's data frame for analysis
m1mf - model.frame(m1)

## If you have rockchalk 1.5.4, run these:
library(rockchalk)
summarize(m1mf)
## otherwise, run
summary(m1mf)

## First, overview for values of xcat1
newdataMaker(m1, fl = list(xcat1 = levels(m1mf$xcat1)))

## mix and match all combinations of xcat1 and xcat2
newdataMaker(m1, fl = list(xcat1 = levels(m1mf$xcat2), xcat2 =
levels(m1mf$xcat2)))

## Pick some particular values for focus
newdataMaker(m1, fl = list(x2 = c(0.25, 1.0), xcat2 = c(M,D)))

## Generate a newdata frame and predictions in one step
predictOMatic(m1, fl = list(x2 = c(0.25, 1.0), xcat2 = c(M,D)))


predictOMatic(m1, fl = list(x2 = c(0.25, 1.0), xcat2 = c(M,D)),
interval=conf)

newdf - predictOMatic(m1, fl = list(x2 = c(0.25, 1.0), xcat2 =
c(M,D), x1=plotSeq(dat$x1)))

plot(y ~ x1, data= model.frame(m1))
by(newdf, list(newdf$x2, newdf$xcat2), function(x) {lines(x$x1, x$fit)})

newdataMaker(m1, fl = list(x2 = c(-1,0, 1), xcat2 = c(M,D)))

predictOMatic(m1, fl = list(x2 = range(dat$x2), xcat2 = c(M,D)))

newdf - predictOMatic(m1, fl = list(x2 = quantile(dat$x2), xcat2 = c(M,D)))
plot(y ~ x2 , data=model.frame(m1))

lines(y ~ x2,  newdf)


predictOMatic(m1, fl = list(x2 = c(50, 60), xcat2 = c(M,D)),
interval=conf)

## just gets the new data
nd - newdataMaker(m1, fl = list(x2 = c(50, 60), xcat2 = c(M,D)))

pr - predictOMatic(m1, fl = list(x2 = c(50, 60), xcat2 = c(M,D)),
interval=conf)

##
m2 - glm(xcat2 ~ x1 + x2 + x3 + xcat1, data=dat, family=binomial(logit))
summary(m2)
m2mf - model.frame(m2)
summarize(m2mf)

predict(m2, newdata=data.frame(xcat1=c(M,F), x1=c(1),
x2=mean(m2mf$x2), x3=mean(m2mf$x3)), se.fit=TRUE)

predictOMatic(m2, fl = list(x2 = c(-1, 1), xcat1 = c(M,F)),
interval=conf, type=response)



Now, the example session ###

 newdataMaker - function (model = NULL, fl = NULL){
+ mf - model.frame(model)
+ mterms - terms(model)
+ mfnames - colnames(mf)
+ modelcv - centralValues(mf)
+ if (sum(!names(fl) %in% mfnames)  0) stop(cat(c(Error. The
focus

Re: [R] Solve an ordinary or generalized eigenvalue problem in R?

[Paul Johnson]

On Fri, Apr 20, 2012 at 2:18 PM, Jonathan Greenberg j...@illinois.edu wrote:
 Ok, I figured out a solution and I'd like to get some feedback on this from
 the R-helpers as to how I could modify the following to be package
 friendly -- the main thing I'm worried about is how to dynamically set the


There's documentation on this in the Writing R Extensions manual,
See 1.2.1 Using Makevars and Section 6.7 Numerical analysis
subrouties. They recommend looking at the package fastICA. It appears
to me you will use some platform neutral Makefile variables.

If you give a full working example from windows--your function, data
to use it--I'll see what I can do on the Linux side.  I've never had
need to worry about this question before, but it is about time I
learned.

Before you think about packaging something like this, I'd say step one
is to clean up your code so it is easier to read.  R CMD check won't
let you get bye with T or F in place of TRUE and FALSE, and you seem
to mix - and = in your code, which makes it difficult to read.

pj


 dyn.load statement below correctly (obviously, its hard coded to my
 particular install, and would only work with windows since I'm using a
 .dll):

 Rdggev - function(JOBVL=F,JOBVR=T,A,B)
 {
 # R implementation of the DGGEV LAPACK function (with generalized
 eigenvalue computation)
 # See http://www.netlib.org/lapack/double/dggev.f
  # coded by Jonathan A. Greenberg j...@illinois.edu
  dyn.load(C:\\Program Files\\R\\R-2.14.2\\bin\\x64\\Rlapack.dll)

 if(JOBVL)
 {
 JOBVL=V
 } else
 {
 JOBVL=N
 }
  if(JOBVR)
 {
 JOBVR=V
 } else
 {
 JOBVR=N
 }
  # Note, no error checking is performed.
  # Input parameters
 N=dim(A)[[1]]
 LDA=N
 LDB=N
 LDVL=N
 LDVR=N
 LWORK=as.integer(max(1,8*N))
  Rdggev_out - .Fortran(dggev, JOBVL, JOBVR, N, A, LDA, B, LDB,
 double(N), double(N), double(N),
 array(data=0,dim=c(LDVL,N)), LDVL, array(data=0,dim=c(LDVR,N)), LDVR,
 double(max(1,LWORK)), LWORK, integer(1))

 names(Rdggev_out)=c(JOBVL,JOBVR,N,A,LDA,B,LDB,ALPHAR,ALPHAI,
 BETA,VL,LDVL,VR,LDVR,WORK,LWORK,INFO)

 Rdggev_out$GENEIGENVALUES=(Rdggev_out$ALPHAR+Rdggev_out$ALPHAI)/Rdggev_out$BETA
  return(Rdggev_out)
 }



 On Fri, Apr 20, 2012 at 12:01 PM, Jonathan Greenberg j...@illinois.eduwrote:

 Incidentally, if you want to test this out:

  A
          [,1]     [,2]     [,3]
 [1,] 1457.738 1053.181 1256.953
 [2,] 1053.181 1213.728 1302.838
 [3,] 1256.953 1302.838 1428.269

  B
          [,1]     [,2]     [,3]
 [1,] 4806.033 1767.480 2622.744
 [2,] 1767.480 3353.603 3259.680
 [3,] 2622.744 3259.680 3476.790

 I THINK this is correct for the other parameters:
         JOBVL=N
         JOBVR=N
 N=dim(A)[[1]]
  LDA=N
 LDB=N
 LDVL=N
  LDVR=N
 LWORK=max(1,8*N)

 .Call(La_dgeev, JOBVL, JOBVR, N, A, LDA, B, LDB, ALPHAR=NA, ALPHAI=NA,
 BETA=NA,
  VL=NA, LDVL, VR=NA, LDVR, WORK=NA, LWORK, INFO=NA, PACKAGE=base)
  LAPACK's documentation is:
 http://www.netlib.org/lapack/double/dggev.f

 --j

 On Fri, Apr 20, 2012 at 11:58 AM, Jonathan Greenberg 
 j...@illinois.eduwrote:

 So I am a complete noob at doing this type of linking.  When I write this
 statement (all the input values are assigned, but I have to confess I don't
 know what to do with the output variables -- in this call they are all
 assigned NA):

 .Call(La_dgeev, JOBVL, JOBVR, N, A, LDA, B, LDB, ALPHAR=NA, ALPHAI=NA,
 BETA=NA,
 + VL=NA, LDVL, VR=NA, LDVR, WORK=NA, LWORK, INFO=NA, PACKAGE=base)

 I get:

 Error in .Call(La_dgeev, JOBVL, JOBVR, N, A, LDA, B, LDB, ALPHAR = NA,
  :
   C symbol name La_dgeev not in DLL for package base

 I'm running this on Windows 7 x64 version of R 2.14.2.  The
 R_ext/Lapack.h file states:

 /* DGEEV - compute for an N-by-N real nonsymmetric matrix A, the */
 /* eigenvalues and, optionally, the left and/or right eigenvectors */
 La_extern void
 F77_NAME(dgeev)(const char* jobvl, const char* jobvr,
 const int* n, double* a, const int* lda,
 double* wr, double* wi, double* vl, const int* ldvl,
  double* vr, const int* ldvr,
 double* work, const int* lwork, int* info);

 Any ideas?

 --j



 On Fri, Apr 20, 2012 at 1:37 AM, Berend Hasselman b...@xs4all.nl wrote:


 On 19-04-2012, at 20:50, Jonathan Greenberg wrote:

  Folks:
 
  I'm trying to port some code from python over to R, and I'm running
 into a
  wall finding R code that can solve a generalized eigenvalue problem
  following this function model:
 
 
 http://docs.scipy.org/doc/scipy/reference/generated/scipy.linalg.eig.html
 
  Any ideas?  I don't want to call python from within R for various
 reasons,
  I'd prefer a native R solution if one exists.  Thanks!

 An old thread is here:
 http://tolstoy.newcastle.edu.au/R/help/05/06/6802.html

 R does provide eigen().
 R has the Lapack routines dggev and dgges in its library.
 You'd have to write your own .Fortran interface to those routines.

 Berend




 --
 Jonathan A. Greenberg, PhD
 Assistant Professor
 Department of Geography and Geographic Information Science
 University of Illinois at Urbana-Champaign
 607 South 

[R] what folder to run write_PACKAGES in?

[Paul Johnson]

I set up a local repo for testing packages. My packages are not
showing up from the repository when viewed by Linux clients. I suspect
this is a web administrator/firewall issue, but it could be I created
the repo wrongly.  I am supposed to run write_PACKAGES separately in
each R-version folder. Right?

Maybe other novices can use these scripts, if they are not wrong :)

Here's the file structure. On the file space that the Web server can
see, I create a folder /tools/kran and directories

 bin
macosx
 leopard
  contrib
2.13
2.14
2.15

windows
 contrib
   2.13
   2.14
   2.15
src
contrib
   2.13
   2.14
   2.15

That's created by this:
#
create_repo_tree - function(local.repos, rversions){

folders  -  c(/bin/windows/contrib,
/bin/macosx/leopard/contrib, /src/contrib)
for(dir in folders){
dirs - paste(local.repos, dir, /, rversions, sep='')
lapply(dirs, dir.create, recursive = TRUE, showWarnings = TRUE)
}
}

create_repo_tree(/tools/kran, c(2.13, 2.14, 2.15))
###

My CRAN mirror is in a sister folder /tools/cran and that works
properly to be served at the address http://rweb.quant.ku.edu/cran.

I want our local testing thing to show at similar
http://rweb.quant.ku.edu/kran.  Supposing the Apache web server magic
is done, I *believe* the following should work.

I dropped packages in the right version folders, and I wrote a script
that goes separately to each version number folder and runs
write_PACKAGES.

### Researchers can upload
### packages into the approrpriate folder.
### Administratively, we schedule this run run every night
write_PACKAGES_wrapper - function(local.repos) {
require(tools)

rversions -  dir(path = paste(local.repos,
/bin/windows/contrib, sep=), full.names = TRUE)
for (i in rversions) write_PACKAGES(dir = i, subdirs=TRUE,
type=win.binary)

#repeat
rversions -  dir(path = paste(local.repos,
/bin/macosx/leopard/contrib, sep=), full.names = TRUE)
for (i in rversions) write_PACKAGES(dir = i, subdirs=TRUE,
type=mac.binary)

rversions -  dir(path = paste(local.repos, /src/contrib,
sep=), full.names = TRUE)
for (i in rversions) write_PACKAGES(dir = i, subdirs=TRUE, type=source)
}


write_PACKAGES_wrapper(/tools/kran)

#

Right?

After running that, I do see the PACKAGES files appear under the
version number directories.

However, from the linux clients I see this:

 install.packages(rockchalk, repos=http://rweb.quant.ku.edu/kran;)
Installing package(s) into ‘/home/pauljohn/R/x86_64-pc-linux-gnu-library/2.15’
(as ‘lib’ is unspecified)
Warning: unable to access index for repository
http://rweb.quant.ku.edu/kran/src/contrib
Warning message:
package ‘rockchalk’ is not available (for R version 2.15.0)

The Web administrator here suggests I've done the write_PACKAGES
incorrectly because there is no PACKAGES file in
/tools/kran/src/contrib.  But I do have PACKAGES files in the
subfolders 2.15.


However, on a windows system, it does work.

 install.packages(rockchalk, repos=http://rweb.quant.ku.edu/kran;)
trying URL 
'http://rweb.quant.ku.edu/kran/bin/windows/contrib/2.15/rockchalk_1.5.5.06.zip'
Content type 'application/zip' length 486682 bytes (475 Kb)
opened URL
downloaded 475 Kb

package ‘rockchalk’ successfully unpacked and MD5 sums checked

The downloaded binary packages are in
C:\Users\pauljohn32\AppData\Local\Temp\Rtmpq0m3Id\downloaded_packages

The Web admin folks say to me, if we did it wrong, nothing would
work. Some does, so it is your fault.

-- 
Paul E. Johnson
Professor, Political Science    Assoc. Director
1541 Lilac Lane, Room 504     Center for Research Methods
University of Kansas               University of Kansas
http://pj.freefaculty.org            http://quant.ku.edu

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

Re: [R] glmmADMB

[Paul Johnson]

On Tue, May 8, 2012 at 5:16 PM, rbuxton moy...@hotmail.com wrote:
 http://r.789695.n4.nabble.com/file/n4618871/Data_for_list_serve.csv
 Data_for_list_serve.csv

 Here is my data, hope this helps.


 The LESP CHUCKLE , FTSP FLIGHT, and ANMU CHIRRUP are the dependent
 variables, I want to run one model for each.


I've not succeeded in replicating your work, the data set has
something funny perhaps.

I saved  your file data.csv and ran this:

callsna - read.table(data.csv, sep=,, header=T)

library(glmmADMB)

mod - glmmadmb(LESP.CHUCKLE ~ Years.Erad + IS + Ref + Dist.Buldir +
Food+Moon+Wind.Speed + (1|SITE/ISLAND), data=callsna,
zeroInflation=TRUE, family=nbinom)

I don't get as far as you do.

 mod - glmmadmb(LESP.CHUCKLE ~ Years.Erad + IS + Ref + Dist.Buldir + 
 Food+Moon+Wind.Speed + (1|SITE/ISLAND), data=callsna, zeroInflation=TRUE, 
 family=nbinom)
Error in II[, ii] = II[, ii] + REmat$codes[[i]] :
  number of items to replace is not a multiple of replacement length
In addition: Warning messages:
1: In glmmadmb(LESP.CHUCKLE ~ Years.Erad + IS + Ref + Dist.Buldir +  :
  NAs removed in constructing fixed-effect model frame: you should
probably remove them manually, e.g. with na.omit()
2: In II[, ii] + REmat$codes[[i]] :
  longer object length is not a multiple of shorter object length
 callsna - na.omit(callsna)
 mod - glmmadmb(LESP.CHUCKLE ~ Years.Erad + IS + Ref + Dist.Buldir + 
 Food+Moon+Wind.Speed + (1|SITE/ISLAND), data=callsna, zeroInflation=TRUE, 
 family=nbinom)
Error in `contrasts-`(`*tmp*`, value = contr.funs[1 + isOF[nn]]) :
  contrasts can be applied only to factors with 2 or more levels


I wondered if the function you are using is supposed to be able to do
that.  I looked back at the web page for glmmADMB (note
capitalization), http://glmmadmb.r-forge.r-project.org/ and it says:

Zero-inflation (currently only as a single constant term across all groups)

But I'll try again if you post your R code that works with that CSV
file you posed, I may try again.

pj

 So, again the desired model is:

 mod - glmmadmb(LESP.CHUCKLE~
 Years.Erad+IS+Ref+Dist.Buldir+Food+Moon+Wind.Speed+(1|SITE/ISLAND),
 data=callsna,
 zeroInflation=TRUE, family=nbinom)

 cheers!
 Rachel


 --
 View this message in context: 
 http://r.789695.n4.nabble.com/glmmADMB-tp4616701p4618871.html
 Sent from the R help mailing list archive at Nabble.com.

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



-- 
Paul E. Johnson
Professor, Political Science    Assoc. Director
1541 Lilac Lane, Room 504     Center for Research Methods
University of Kansas               University of Kansas
http://pj.freefaculty.org            http://quant.ku.edu

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

Re: [R] low R square value from ANCOVA model

[Paul Johnson]

On Tue, May 8, 2012 at 3:45 PM, array chip arrayprof...@yahoo.com wrote:
 Thanks again Peter. What about the argument that because low R square (e.g. 
 R^2=0.2) indicated the model variance was not sufficiently explained by the 
 factors in the model, there might be additional factors that should be 
 identified and included in the model. And If these additional factors were 
 indeed included, it might change the significance for the factor of interest 
 that previously showed significant coefficient. In other word, if R square is 
 low, the significant coefficient observed is not trustworthy.

 What's your opinion on this argument?

I think that argument is silly. I'm sorry if that is too blunt. Its
just plain superficial.
 It reflects a poor understanding of what the linear model is all
about. If you have
other variables that might belong in the model, run them and test.
The R-square,
either low or high, does not have anything direct to say about whether
those other
variables exist.

Here's my authority.

Arthur Goldberger (A Course in Econometrics, 1991, p.177)
“Nothing in the CR (Classical Regression) model requires that R2 be high. Hence,
a high R2 is not evidence in favor of the model, and a low R2 is not evidence
against it.”

I found that reference in Anders Skrondal and  Sophia Rabe-Hesketh,
Generalized Latend Variable Modeling: Multilevel, Longitudinal,
and Structural Equation Models, Boca Raton, FL: Chapman and Hall/CRC, 2004.

From Section 8.5.2:

Furthermore, how badly the baseline model fits the data depends greatly
on the magnitude of the parameters of the true model. For instance, consider
estimating a simple parallel measurement model. If the true model is a
congeneric measurement model (with considerable variation in factor loadings
and measurement error variances between items), the fit index could be high
simply because the null model fits very poorly, i.e. because the
reliabilities of
the items are high. However, if the true model is a parallel measurement model
with low reliabilities the fit index could be low although we are estimating the
correct model. Similarly, estimating a simple linear regression model can yield
a high R2 if the relationship is actually quadratic with a considerable linear
trend and a low R2 when the model is true but with a small slope (relative to
the overall variance).

For a detailed argument/explanation of the argument that the R-square is not
a way to decide if a model is good or bad see

King, Gary. (1986). How Not to Lie with Statistics: Avoiding Common Mistakes in
Quantitative Political Science. American Journal of Political Science,
30(3), 666–687. doi:10.2307/2111095

pj
-- 
Paul E. Johnson
Professor, Political Science    Assoc. Director
1541 Lilac Lane, Room 504     Center for Research Methods
University of Kansas               University of Kansas
http://pj.freefaculty.org            http://quant.ku.edu

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

Re: [R] glmmADMB

[Paul Johnson]

On Tue, May 8, 2012 at 9:06 PM, rbuxton moy...@hotmail.com wrote:

 Update!

 I changed the site categories.  I noticed that I had coded them as
 North,
 South, East, West on different islands, which may have caused confusion in
 the model.

 [...]



 mod - glmmadmb(LESP.CHUCKLE~ 1+(1|ISLAND), data=callsna,
 zeroInflation=TRUE, family=nbinom)

 Any thoughts?

 YES.

POST THE ENTIRE R PROGRAM YOU ARE RUNNING as well as a link to the EXACT
data set that is causing the problem.

pj




 Thanks so so much!

 Rachel Buxton





-- 
Paul E. Johnson
Professor, Political ScienceAssoc. Director
1541 Lilac Lane, Room 504 Center for Research Methods
University of Kansas   University of Kansas
http://pj.freefaculty.orghttp://quant.ku.edu

[[alternative HTML version deleted]]

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

[R] Windows Task Scheduler and R updates. Need basic tips

[Paul Johnson]

This is a basic Windows system administrator problem, asked by a Linux
guy who is helping out in a Windows lab.

I want to keep R packages up to date on MS Windows 7 with a job in the
Task Scheduler.  I have an R program that I can run (as
administrator) that updates the existing packages and then installs
all the new ones.

I do not understand how to run that in a dependable way in the scheduler.

If I put the update script R-update.R in, for example, in

C:\Program Files\R\R-update.R

Then what?  Do I need a CMD batch script to run the R script?

I can't tell where Windows wants to write the standard output and
error for my R job.

And while I'm asking, does Windows care if I run

R CMD BATCH C:\Program Files\R\R-update.R

or

R --vanilla -f C:\Program Files\R\R-update.R

pj
-- 
Paul E. Johnson
Professor, Political Science    Assoc. Director
1541 Lilac Lane, Room 504     Center for Research Methods
University of Kansas               University of Kansas
http://pj.freefaculty.org            http://quant.ku.edu

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

Re: [R] Names of Greek letters stored as character strings; plotmath.

[Paul Johnson]

On Sat, May 19, 2012 at 11:07 AM, William Dunlap wdun...@tibco.com wrote:
 parse(text=paste(...)) works in simple cases but not in others.  The
 fortune about it is there because it is tempting to use but if you bury it
 in a general purpose function it will cause problems when people
 start using nonstandard names for variables.  bquote(), substitute(),
 call(), and relatives work in all cases.  E.g.,

   par(mfrow=c(2,1))
   power - gamma ; x - Waist ; y - Weight # valid R variable names
   plot(0, main=bquote(.(as.name(x))^.(as.name(power))/.(as.name(y
   plot(0, main=parse(text=paste0(x, ^, power, /, y))) # same as previous
  
   power - gamma ; x - Waist Size (cm) ; y - Weight (kg) # invalid R 
 names
   plot(0, main=bquote(.(as.name(x))^.(as.name(power))/.(as.name(y
   plot(0, main=parse(text=paste0(x, ^, power, /, y))) # whoops
  Error in parse(text = paste0(x, ^, power, /, y)) :
    text:1:7: unexpected symbol
  1: Waist Size
           ^

 Now you might say that serves me right for using weird variable names,
 but some of us use R as a back end to a GUI system (one not designed
 around R) and don't want to inflict on users R's rules for names when
 we do not have to.

 Bill Dunlap
 Spotfire, TIBCO Software
 wdunlap tibco.com

Bill's point is on the money here. We should aim to teach one way that
works always. But finding that one way is the hard part for me.

Lately, I'm bothered that R (or computers in general?) allows too many
ways to do the same thing that work SOME of the time. Without a very
deep understanding of the terminology, users are bewildered.

Go read ?plotmath.  Do we see as.name?   NO.  Why does the idiom
expression(paste()) work as well?   (Don't answer, that's a rhetorical
question).

There are too many ways to do the same thing. Or, well, too many of us
know one way that works sometime and don't find out it doesn't always
work until, as Bill says, it is buried in the bottom of some big
function that behaves erratically.

Gabor notes this works (sometimes):

plot(0, xlab = parse(text = xNm))

Bert observes this works (sometimes)

xnm - quote(gamma)
## makes xnm the name gamma not the string gamma
plot(0,xlab = bquote( .(xnm))

Bill observes this works:
xnm - gamma
plot(0,xlab = bquote(.(as.name(xnm

In line with that, Bill suggests:

power - gamma ; x - Waist Size (cm) ; y - Weight (kg) # invalid R names
plot(0, main=bquote(.(as.name(x))^.(as.name(power))/.(as.name(y

It appears to me that 2/3 of the as.name usages are not needed, at
least not in R 2.15 on Debian Linux.

This works just as well

plot(0, main=bquote(.(x)^.(as.name(power))/.(y)))


Is there any important reason to wrap x and y in this command with as.name?

It does work, maybe I be in the general habit of wrapping as.name
around variables in bquotes?  I'm more confused.

The good new for me is that I cannot find any replacements for the
usage of as.name in that particular expression. So there is just one
way.

Oh, wait, I spoke too soon. as.symbol is identical to as.name.

plot(0, main=bquote(.(x)^.(as.symbol(power))/.(y)))

And then, logically, this works:

plot(0, main=bquote(.(x)^.(as.vector(power, symbol))/.(y)))

I personally think as.symbol is more likely to be understood by
students, I may stick to that. So I think the most succinct, best
approach is


plot(0, main=bquote(.(x)^.(as.symbol(power))/.(y)))

or perhaps the most consistent-looking strategy is:

power - as.symbol(gamma) ; x - Waist Size (cm) ; y - Weight (kg)
plot(0, main=bquote(.(x)^.(power)/.(y)))

Speaking of commands that have identical results and exist separately
for no apparently good reason:

library(help=rockchalk)

help(package=rockchalk)

pj


 -Original Message-
 From: r-help-boun...@r-project.org [mailto:r-help-boun...@r-project.org] On 
 Behalf
 Of Bert Gunter
 Sent: Saturday, May 19, 2012 7:24 AM
 To: Gabor Grothendieck
 Cc: r-help
 Subject: Re: [R] Names of Greek letters stored as character strings; 
 plotmath.

 ... and here is another incantation that may be  informative.

 xnm- as.name(gamma')  ## This does the parsing
 plot(0, xlab =bquote(.(xnm))

 The initial puzzle is that if you just set
 xnm - gamma

 bquote will insert the string gamma rather than the symbol. After
 all, that's what plotmath sees for xnm. So the key is telling plotmath
 that it's a symbol, not a string. This can either be done before, as
 above, or inline, as you and Gabor showed. Unsurprisingly. this also
 does it, since as.name() is doing the parsing:

 xnm - gamma
  plot(0,xlab=bquote(.(as.name(xnm

 AND we are adhering to Thomas's dictum: bquote is a wrapper for
 substitute(), which is what he recommends as the preferable
 alternative to eval(parse(...)) . But, heck -- all such software
 principles are just guidelines. Whatever works (robustly).

 HTH.

 Cheers,
 Bert

 On Sat, May 19, 2012 at 3:17 AM, Gabor Grothendieck
 ggrothendi...@gmail.com wrote:
  On Sat, May 19, 2012 at 1:18 AM, Rolf Turner rolf.tur...@xtra.co.nz 
 

[R] Build error on RedHat EL 5.5: byte-compiling package 'compiler'

[Paul Johnson]

Greetings R-help!

In case anybody has worked on an old Redhat system lately, can you
remember what is the cause of this problem?  I don't think this is a
fatal problem, because I can get around it by uninstalling all of the
R RPM packages and re-running the build. But if R is installed, the
build fails while byte-compiling as seen below.  A *guessed* removing
the installed R stuff might make a difference, and it did!

The build directory is /tmp/BUILD/R-2.15.0/ but the compile process
jumps out and tries to grab files in /usr/lib64...


make[6]: Entering directory `/tmp/BUILD/R-2.15.0/src/library/tools/src'
make[6]: `Makedeps' is up to date.
make[6]: Leaving directory `/tmp/BUILD/R-2.15.0/src/library/tools/src'
make[6]: Entering directory `/tmp/BUILD/R-2.15.0/src/library/tools/src'
mkdir -p -- ../../../../library/tools/libs
make[6]: Leaving directory `/tmp/BUILD/R-2.15.0/src/library/tools/src'
make[5]: Leaving directory `/tmp/BUILD/R-2.15.0/src/library/tools/src'
make[4]: Leaving directory `/tmp/BUILD/R-2.15.0/src/library/tools'
make[3]: Leaving directory `/tmp/BUILD/R-2.15.0/src/library/tools'
make[3]: Entering directory `/tmp/BUILD/R-2.15.0/src/library/compiler'
building package 'compiler'
mkdir -p -- ../../../library/compiler
make[4]: Entering directory `/tmp/BUILD/R-2.15.0/src/library/compiler'
mkdir -p -- ../../../library/compiler/R
mkdir -p -- ../../../library/compiler/po
make[4]: Leaving directory `/tmp/BUILD/R-2.15.0/src/library/compiler'
make[4]: Entering directory `/tmp/BUILD/R-2.15.0/src/library/compiler'
byte-compiling package 'compiler'
Warning in file(datafile, wb) :
  cannot open file '/usr/lib64/R/library/compiler/R/compiler.rdb':
Permission denied
Error in file(datafile, wb) : cannot open the connection
Calls: Anonymous - code2LazyLoadDB - makeLazyLoadDB - close - file
Execution halted
make[4]: *** [../../../library/compiler/R/compiler.rdb] Error 1
make[4]: Leaving directory `/tmp/BUILD/R-2.15.0/src/library/compiler'
make[3]: *** [all] Error 2
make[3]: Leaving directory `/tmp/BUILD/R-2.15.0/src/library/compiler'
make[2]: *** [R] Error 1
make[2]: Leaving directory `/tmp/BUILD/R-2.15.0/src/library'
make[1]: *** [R] Error 1
make[1]: Leaving directory `/tmp/BUILD/R-2.15.0/src'
make: *** [R] Error 1
error: Bad exit status from /var/tmp/rpm-tmp.5874 (%build)


The command line that initiates the RPM build is

rpmbuild -ba R-2.15-pj.spec --define dist=kurocks  Rbuild.txt 21 

and I uploaded the full script of the build process here:

http://pj.freefaculty.org/R/Rbuild.txt

I searched the archive and find a 2006 mention of the same problem,
and the solution was that the user had R_LIBS defined in the
environment.  I don't have anything like that in the output of env,
but obviously, I'm getting something that confuses the build from the
R envionment.  I am just a little curious if somebody can explain
what's going wrong.

pj
-- 
Paul E. Johnson
Professor, Political Science    Assoc. Director
1541 Lilac Lane, Room 504     Center for Research Methods
University of Kansas               University of Kansas
http://pj.freefaculty.org            http://quant.ku.edu

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

Re: [R] R Memory Issues

[Paul Johnson]

Dear Emiliano:

When they say to read the posting guide, mostly they mean read the
posting guide. But I'll tell you the short version.

1. Include a full runable R program that causes the trouble you are
concerned about.  Include the data or a link to the data, usually the
smallest possible example is what they want.  They don't want 1000
lines of your dissertation project, they want 10 lines needed to
produce the problem you are concerned about.

The point here is this: Don't make people guess about what commands
you ran or what your data actually was. You are going to get the
attention of these folks one time, and you waste it by not reading the
guide and not giving the full details.

2. Include the output from sessionInfo() whenever you ask a question
of this sort.


 sessionInfo()
R version 2.15.0 (2012-03-30)
Platform: x86_64-pc-linux-gnu (64-bit)

locale:
 [1] LC_CTYPE=en_US.UTF-8   LC_NUMERIC=C
 [3] LC_TIME=en_US.UTF-8LC_COLLATE=en_US.UTF-8
 [5] LC_MONETARY=en_US.UTF-8LC_MESSAGES=en_US.UTF-8
 [7] LC_PAPER=C LC_NAME=C
 [9] LC_ADDRESS=C   LC_TELEPHONE=C
[11] LC_MEASUREMENT=en_US.UTF-8 LC_IDENTIFICATION=C

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


I can tell  you my best guess about what is wrong: I suspect you have
a corrupted R install. If you had given us the full working code, I
could have tested that theory. But, alas, I can't.

Why do I think so? I've taught a course this term with 45 students and
about 1 time per week, a student would turn up with that can't
allocate vector... error you see.  On Windows, sometimes it seems the
problem is due to installing R as an administrator and then trying to
update some packages as a non-administrator.  In one really
frustrating case, student has installed car both as admin and as the
user, and the one that was at the front of the search path was
damaged, but we kept removing and re-installing the other one and
nothing was fixed.  Until I noticed there were 2

pj

On Tue, May 22, 2012 at 11:40 AM, Emiliano Zapata ezapata...@gmail.com wrote:
 As a continuation to my original question, here is the massage that I get:

 Error in glm.fit(x = structure(c(1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1,  :
  cannot allocate memory block of size 2.1 Gb

 The model glm.fit is a logistic type (in the family of GLM) model. Maybe
 this is not enough information; again!, but some feedback will be
 appreciated. To me the issues appears to be associated with manipulation of
 large dataset. Howeverl the algorithm runs fine in Unix; but not in Windows
 (64 bits windows 7).

 EZ

 On Sun, May 20, 2012 at 4:09 PM, Emiliano Zapata ezapata...@gmail.comwrote:

 Already then, thank you everyone. This information was extremly useful,
 and I'll do a better job on the web next time.

 On Sun, May 20, 2012 at 2:10 PM, Prof Brian Ripley 
 rip...@stats.ox.ac.ukwrote:

 On 20/05/2012 18:42, jim holtman wrote:

 At the point in time that you get the error message, how big are the
 objects that you have in memory?  What does 'memory.size()' show as
 being used?  What does 'memory.limit()' show?  Have you tried using
 'gc()' periodically to do some garbage collection?  It might be that
 you memory is fragmented.  You need to supply some additional
 information.


 Either this is a 32-bit version of R in which case the wrong version is
 being used, or your advice is wrong: there are no credible fragmentation
 issues (and no need to use gc()) on a 64-bit build of R.

 But, we have a posting guide, we require 'at a minimum information', and
 the OP failed to give it to us so we are all guessing, completely
 unnecessarily.



 On Sun, May 20, 2012 at 12:09 PM, Emiliano Zapataezapata...@gmail.com
  wrote:

 -- Forwarded message --
 From: Emiliano Zapataezapata...@gmail.com
 Date: Sun, May 20, 2012 at 12:09 PM
 Subject:
 To: R-help@r-project.org


 Hi,

 I have a 64 bits machine (Windows) with a total of 192GB of physical
 memory
 (RAM), and total of 8 CPU. I wanted to ask how can I make R make use of
 all
 the memory. I recently ran a script requiring approximately 92 GB of
 memory
 to run, and got the massage:

  cannot allocate memory block of size 2.1 Gb



 I read on the web that if you increase the memory you have to reinstall
 R;
 would that be enough. Could I just increase the memory manually.


 Take you for any comments, or links on the web.


 EZ

        [[alternative HTML version deleted]]

 __**
 R-help@r-project.org mailing list
 https://stat.ethz.ch/mailman/**listinfo/r-helphttps://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.






 --
 Brian D. Ripley,                  rip...@stats.ox.ac.uk
 Professor of Applied Statistics,  
 

[R] Ask if an object will respond to a function or method

[Paul Johnson]

In the rockchalk package, I want to provide functions for regression
objects that are "well behaved." If an object responds to the methods
that lm or glm objects can handle, like coef(), nobs(), and summary(),
I want to be able to handle the same thing.

It is more difficult than expected to ask a given fitted model object
"do you respond to these functions: coef(), nobs(), summary()." How
would you do it?

I tried this with the methods() function but learned that all methods
that a class can perform are not listed.  I'll demonstrate with a
regression "zz" that is created by the example in the plm package.
The coef() function succeeds on the zz object, but coef is not listed
in the list of methods that the function can carry out.

> library(plm)
> example(plm)

> class(zz)
[1] "plm""panelmodel"
> methods(class = "plm")
 [1] ercomp  fixef   has.intercept   model.matrix
 [5] pFtest  plmtest plotpmodel.response
 [9] pooltestpredict residuals   summary
[13] vcovBK  vcovDC  vcovG   vcovHC
[17] vcovNW  vcovSCC
see '?methods' for accessing help and source code
> methods(class = "panelmodel")
 [1] deviance  df.residual   fittedhas.intercept index
 [6] nobs  pbgtest   pbsytest  pcdtest   pdim
[11] pdwtest   phtestprint pwartest  pwfdtest
[16] pwtestresiduals terms updatevcov
see '?methods' for accessing help and source code
> coef(zz)
   log(pcap)  log(pc) log(emp)unemp
-0.026149654  0.292006925  0.768159473 -0.005297741

I don't understand why coef(zz) succeeds but coef is not listed as a method.

Right now, I'm contemplating this:

zz1 < - try(coef(zz))
if (inherits(zz1, "try-error")) stop("Your model has no coef method")

This seems like a bad workaround because I have to actually run the
function in order to find out if the function exists. That might be
time consuming for some summary() methods.

pj

-- 
Paul E. Johnson
Professor, Political ScienceDirector
1541 Lilac Lane, Room 504  Center for Research Methods
University of Kansas University of Kansas
http://pj.freefaculty.org  http://crmda.ku.edu

__
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] xvfb? cron job updates R packages, fails on some requiring X11

[Paul Johnson]

In Centos 7 systems, I wrote a script that runs on the cron and I
notice some package updates and installs fail like this:


Error : .onLoad failed in loadNamespace() for 'iplots', details:
  call: .jnew("org/rosuda/iplots/Framework")
  error: java.awt.HeadlessException:
No X11 DISPLAY variable was set, but this program performed an
operation which requires it.
Error: loading failed
Execution halted
ERROR: loading failed
* removing ‘/usr/share/R/library/iplots’
* restoring previous ‘/usr/share/R/library/iplots’


I can log in interactively and run  the same install successfully.

I understand I need something like xvfb to simulate an X11 session,
but I don't understand how to make it work.  Can one of you give me an
idiot's guide on what to do?

pj


-- 
Paul E. Johnson   http://pj.freefaculty.org
Director, Center for Research Methods and Data Analysis http://crmda.ku.edu

To write to me directly, please address me at pauljohn at ku.edu.

__
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 installing packages

[Paul Johnson]

On Wed, Oct 19, 2016 at 10:31 AM, David Winsemius
 wrote:
>
>> On Oct 19, 2016, at 4:54 AM, Kevin E. Thorpe  
>> wrote:
>>
>> Hello.
>>
>> I am posting this on behalf of one of my students who is getting error 
>> messages when installing some packages. I have not seen this before nor have 
>> I been able to replicate it. I'm including the relevant (I think) 
>> information. I get the students to install rms with dependencies. As you can 
>> see, rms does get installed but when the attempt is made to attach it, 
>> ggplot2 cannot be loaded. Thus I tried explicitly installing ggplot2 and you 
>> can see then ensuing errors below. I have included the sessionInfo() at the 
>> end.
>>
>> I hope someone can point me at a solution.
>>

I have just seen this on Win10 for the first time on a student's new
laptop.  Anti-virus software is possible reason for this, but I don't
know.

I did find a very direct solution.

WHen you run the package install, and it says "cannot move from
temporary C:\Users\yourname\AppData\XYZ to user directory
C:\Users\yourname\Documents\R\site-library", just move the folder
manually.  It works every time.

Unfortunately, Windows hides AppData.  In the Windows explorer view
options, tell it not to hide protected files.  Then in explorer
navigate into the user's AppData\Temp\Rxxx folder, you'll see the
downloaded zip files there for Rcpp and such.  Doubleclick the zip
file, right click copy the directory name "Rcpp" and paste it into the
user's R folder, under C:/Users/yourname/Documents/R/3.3/site-library
(or whatever that's called).

We ran into about 2 packages that have this failure to copy from
temporary space, but this old-fashioned copy a folder over method
worked fine. After this, R was able to load Rcpp, RcppEigen.

There is no indication that those R files are in a virus quarantine,
so I can't say for sure the security software is the cause.  At first,
I thought the problem was the usual one that we have been aware of for
some time--Windows thinks those files are in use and will not replace
them.

 We are trying to remove a layer of virus protection programs
installed by Dell to see if this happens later. If I hear back from
that student, I'll let you know. Another layer in this story is that
her 30 day free trials were expired, and I have feeling this means
that not only is McAfee still installed, but it is also running but
refusing to let you interact with it.

pj

-- 
Paul E. Johnson   http://pj.freefaculty.org
Director, Center for Research Methods and Data Analysis http://crmda.ku.edu

I only use this account for email list memberships. To write directly,
address me at pauljohn at ku.edu.

__
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] Interesting quirk with fractions and rounding

[Paul Johnson]

Hello, R friends

My student unearthed this quirk that might interest you.

I wondered if this might be a bug in the R interpreter. If not a bug,
it certainly stands as a good example of the dangers of floating point
numbers in computing.

What do you think?

> 100*(23/40)
[1] 57.5
> (100*23)/40
[1] 57.5
> round(100*(23/40))
[1] 57
> round((100*23)/40)
[1] 58

The result in the 2 rounds should be the same, I think.  Clearly some
digital number devil is at work. I *guess* that when you put in whole
numbers and group them like this (100*23), the interpreter does
integer math, but if you group (23/40), you force a fractional
division and a floating point number. The results from the first 2
calculations are not actually 57.5, they just appear that way.

Before you close the books, look at this:

> aa <- 100*(23/40)
> bb <- (100*23)/40
> all.equal(aa,bb)
[1] TRUE
> round(aa)
[1] 57
> round(bb)
[1] 58

I'm putting this one in my collection of "difficult to understand"
numerical calculations.

If you have seen this before, I'm sorry to waste your time.

pj
-- 
Paul E. Johnson   http://pj.freefaculty.org
Director, Center for Research Methods and Data Analysis http://crmda.ku.edu

To write to me directly, please address me at pauljohn at ku.edu.

__
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] Prediction plots

[Paul Johnson]

I have done this a lot. Would you mind installing my pkg rockchalk and then
run example(plotSlope) and example(plotCurve)? If the output is close to
what you want, you can adjust my code. The vignette explains.

1. Create newdata object
2. Run that through predict
3. Make plot

None of this is rocket science, but it does help students here.

PJ

On Apr 18, 2017 1:17 PM, "Santiago Bueno"  wrote:

> Thanks Boris, the following is an extract of my data. I have developed
> biomass models using codes like:
>
> start <- coef (lm(log(Btot)~I(log(dbh**2*haut)),data=dat[dat$Btot>0,]))
>
> start[1] <- exp(start[1])
>
> names(start) <- c("a","b")
>
> M1 <- nls(Btot~a*(dbh**2*haut)**b,data=dat,start=start,weights=
> 1/dat$dbh**4)
>
>
> start <- coef(lm(log(Btot)~I(log(dbh))+I(log(haut)),data=dat[dat$
> Btot>0,]))
>
> start[1] <- exp(start[1])
>
> names(start) <- c("a","b1","b2")
>
> M2 <- nls(Btot~a*dbh**b1*haut**b2,data=dat,start=start,weights=
> 1/dat$dbh**4)
>
>
> Tree No dbh haut Btot
> 1 35.00 18.90 0.535
> 2 25.00 16.60 0.248
> 3 23.00 19.50 0.228
> 4 13.50 15.60 0.080
> 5 20.00 18.80 0.172
> 6 23.00 17.40 0.190
> 7 29.00 19.90 0.559
> 8 17.60 18.20 0.117
> 9 31.00 25.30 0.645
> 10 26.00 23.50 0.394
> 11 13.00 13.00 0.038
> 12 32.00 20.70 0.443
> It is my interest to get prediction plots for the models. I have tried to
> use the following codes with no success: Let m be one of the fitted models
> with dbh as the only entry. To construct a plot of the predictions made by
> this model I have tried:
> with(dat,plot(dbh,Btot,xlab="Dbh(cm)",ylab="Biomass (t)"))
> D <- seq(par("usr")[1],par("usr")[2],length=200)
> lines(D,predict(m,newdata=data.frame(dbh=D)),col="red")
> For a model m that has dbh and height as entries, I have tried to get its
> predictions as follows:
> D <- seq(0,180,length=20)
> H <- seq(0,61,length=20)
> B <- matrix(predict(m,newdata=expand.grid(dbh=D,height=H)),length(D))
>
> Can someone provide help please!!!
>
>
> Best regards,
>
> Santiago Bueno
>
> [[alternative HTML version deleted]]
>
> __
> R-help@r-project.org mailing list -- To UNSUBSCRIBE and more, see
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide http://www.R-project.org/
> posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
>

[[alternative HTML version deleted]]

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

Re: [R] Interesting quirk with fractions and rounding

[Paul Johnson]

We all agree it is a problem with digital computing, not unique to R. I
don't think that is the right place to stop.

What to do? The round example arose in a real funded project where 2 R
programs differed in results and cause was  that one person got 57 and
another got 58. The explanation was found, but its less clear how to
prevent similar in future. Guidelines, anyone?

So far, these are my guidelines.

1. Insert L on numbers to signal that you really mean INTEGER. In R,
forgetting the L in a single number will usually promote whole calculation
to floats.
2. S3 variables are called 'numeric' if they are integer or double storage.
So avoid "is.numeric" and prefer "is.double".
3. == is a total fail on floats
4. Run print with digits=20 so we can see the less rounded number. Perhaps
start sessions with "options(digits=20)"
5. all.equal does what it promises, but one must be cautious.

Are there math habits we should follow?

For example, Is it generally true in R that (100*x)/y is more accurate than
100*(x/y), if x > y?   (If that is generally true, couldn't the R
interpreter do it for the user?)

I've seen this problem before. In later editions of the game theory program
Gambit, extraordinary effort was taken to keep values symbolically as
integers as long as possible. Avoid division until the last steps. Same in
Swarm simulations. Gary Polhill wrote an essay about the Ghost in the
Machine along those lines, showing accidents from trusting floats.

I wonder now if all uses of > or < with numeric variables are suspect.

Oh well. If everybody posts their advice, I will write a summary.

Paul Johnson
University of Kansas

On Apr 21, 2017 12:02 AM, "PIKAL Petr" <petr.pi...@precheza.cz> wrote:

> Hi
>
> The problem is that people using Excel or probably other such spreadsheets
> do not encounter this behaviour as Excel silently rounds all your
> calculations and makes approximate comparison without telling it does so.
> Therefore most people usually do not have any knowledge of floating point
> numbers representation.
>
>  Cheers
> Petr
>
> -----Original Message-
> From: R-help [mailto:r-help-boun...@r-project.org] On Behalf Of Paul
> Johnson
> Sent: Thursday, April 20, 2017 11:56 PM
> To: R-help <r-help@r-project.org>
> Subject: [R] Interesting quirk with fractions and rounding
>
> Hello, R friends
>
> My student unearthed this quirk that might interest you.
>
> I wondered if this might be a bug in the R interpreter. If not a bug, it
> certainly stands as a good example of the dangers of floating point numbers
> in computing.
>
> What do you think?
>
> > 100*(23/40)
> [1] 57.5
> > (100*23)/40
> [1] 57.5
> > round(100*(23/40))
> [1] 57
> > round((100*23)/40)
> [1] 58
>
> The result in the 2 rounds should be the same, I think.  Clearly some
> digital number devil is at work. I *guess* that when you put in whole
> numbers and group them like this (100*23), the interpreter does integer
> math, but if you group (23/40), you force a fractional division and a
> floating point number. The results from the first 2 calculations are not
> actually 57.5, they just appear that way.
>
> Before you close the books, look at this:
>
> > aa <- 100*(23/40)
> > bb <- (100*23)/40
> > all.equal(aa,bb)
> [1] TRUE
> > round(aa)
> [1] 57
> > round(bb)
> [1] 58
>
> I'm putting this one in my collection of "difficult to understand"
> numerical calculations.
>
> If you have seen this before, I'm sorry to waste your time.
>
> pj
> --
> Paul E. Johnson   http://pj.freefaculty.org
> Director, Center for Research Methods and Data Analysis
> http://crmda.ku.edu
>
> To write to me directly, please address me at pauljohn at ku.edu.
>
> __
> 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.
>
> 
> Tento e-mail a jakékoliv k němu připojené dokumenty jsou důvěrné a jsou
> určeny pouze jeho adresátům.
> Jestliže jste obdržel(a) tento e-mail omylem, informujte laskavě
> neprodleně jeho odesílatele. Obsah tohoto emailu i s přílohami a jeho kopie
> vymažte ze svého systému.
> Nejste-li zamýšleným adresátem tohoto emailu, nejste oprávněni tento email
> jakkoliv užívat, rozšiřovat, kopírovat či zveřejňovat.
> Odesílatel e-mailu neodpovídá za eventuální škodu způsobenou modifikacemi
> či zpožděním přenosu e-mailu.
>
> V případě, že je tento e-mail součástí obchodního jednání:
> - vyhrazuje si odesílate

Re: [R] Peformance question

[Paul Johnson]

I dont understand your code. But I do have suggestion. Run the functions in
the profiler, maybe differences will point at the enemy.

Know what I mean?

Rprof('check.out')
#run code
Rprof(NULL)
summaryRprof('check.out')

Do that for each method. That may be uninformative.

I wondered if you tried to compile your functions? In some cases it helps
erase differences like this. Norman Matloff has examples like that in Art
of R Programming.

I keep a list of things that are slow, if we can put finger on problem, I
will add to list. I suspect slow here is in runtime object lookup. The
environment ones have info located more quickly by the runtime, I expect.
Also, passing info back and forth from the R runtime system using [ is a
common cause of slow. It is why everybody yells 'vectorize' and 'use
lapply' all the time.  Then again, I'm guessing because I dont understand
your code:)

Good luck,
PJ


On Apr 11, 2017 7:44 PM, "Thomas Mailund"  wrote:

Hi y’all,

I’m working on a book on how to implement functional data structures in R,
and in particular on a chapter on implementing queues. You get get the
current version here https://www.dropbox.com/s/9c2yk3a67p1ypmr/book.pdf?dl=0
and the relevant pages are 50-59. I’ve implemented three versions of the
same idea, implementing a queue using two linked lists. One list contains
the elements you add to the end of a list, the other contains the elements
at the front of the list, and when you try to get an element from a list
and the front-list is empty you move elements from the back-list to the
front. The asymptotic analysis is explained in this figure
https://www.dropbox.com/s/tzi84zmyq16hdx0/queue-
amortized-linear-bound.png?dl=0 and all my implementations do get a linear
time complexity when I evaluate them on a linear number of operations.
However, the two implementations that uses environments seem to be almost
twice as fast as the implementation that gives me a persistent data
structure (see https://www.dropbox.com/s/i9dyab9ordkm0xj/queue-
comparisons.png?dl=0), and I cannot figure out why.

The code below contains the implementation of all three versions of the
queue plus the code I use to measure their performances. I’m sorry it is a
little long, but it is a minimal implementation of all three variants, the
comments just make it look longer than it really is.

Since the three implementations are doing basically the same things, I am a
little stumped about why the performance is so consistently different.

Can anyone shed some light on this, or help me figure out how to explore
this further?

Cheers

Thomas



## Implementations of queues ##

#' Test if a data structure is empty
#' @param x The data structure
#' @return TRUE if x is empty.
#' @export
is_empty <- function(x) UseMethod("is_empty")

#' Add an element to a queue
#' @param x A queue
#' @param elm An element
#' @return an updated queue where the element has been added
#' @export
enqueue <- function(x, elm) UseMethod("enqueue")

#' Get the front element of a queue
#' @param x A queue
#' @return the front element of the queue
#' @export
front <- function(x) UseMethod("front")

#' Remove the front element of a queue
#' @param x The queue
#' @return The updated queue
#' @export
dequeue <- function(x) UseMethod("dequeue")

## Linked lists #

#' Add a head item to a linked list.
#' @param elem  The item to put at the head of the list.
#' @param lst   The list -- it will become the tail of the new list.
#' @return a new linked list.
#' @export
list_cons <- function(elem, lst)
  structure(list(head = elem, tail = lst), class = "linked_list")

list_nil <- list_cons(NA, NULL)

#' @method is_empty linked_list
#' @export
is_empty.linked_list <- function(x) identical(x, list_nil)

#' Create an empty linked list.
#' @return an empty linked list.
#' @export
empty_list <- function() list_nil


#' Get the item at the head of a linked list.
#' @param lst The list
#' @return The element at the head of the list.
#' @export
list_head <- function(lst) lst$head

#' Get the tail of a linked list.
#' @param lst The list
#' @return The tail of the list
#' @export
list_tail <- function(lst) lst$tail

#' Reverse a list
#' @param lst A list
#' @return the reverse of lst
#' @export
list_reverse <- function(lst) {
  acc <- empty_list()
  while (!is_empty(lst)) {
acc <- list_cons(list_head(lst), acc)
lst <- list_tail(lst)
  }
  acc
}


## Environment queues #

queue_environment <- function(front, back) {
  e <- new.env(parent = emptyenv())
  e$front <- front
  e$back <- back
  class(e) <- c("env_queue", "environment")
  e
}

#' Construct an empty closure based queue
#' @return an empty queue
#' @export
empty_env_queue <- function()
  queue_environment(empty_list(), empty_list())

#' @method is_empty env_queue
#' @export
is_empty.env_queue <- function(x)
  is_empty(x$front) && is_empty(x$back)

#' @method enqueue 

Re: [R] Interesting quirk with fractions and rounding

[Paul Johnson]

On Apr 21, 2017 12:01 PM, "JRG" <loesl...@accucom.net> wrote:

A good part of the problem in the specific case you initially presented
is that some non-integer numbers have an exact representation in the
binary floating point arithmetic being used.  Basically, if the
fractional part is of the form 1/2^k for some integer k > 0, there is an
exact representation in the binary floating point scheme.

> options(digits=20)

> (100*23)/40
[1] 57.5

> 100*(23/40)
[1] 57.492895

So the two operations give a slightly different result because the
fractional part of the division of 100*23 by 40 is 0.5.  So the first
operations gives, exactly, 57.5 while the second operation does not
because 23/40 has no exact representation.

Thanks for answering.

This case seemed fun because it was not a contrived example.  We found this
one by comparing masses of report tables from 2 separate programs. It
happened 1 time in about 10,000 calculations.

Guidelines for R coders, though, would be welcome. So far, all I am sure of
is

1 Don't use == for floating point numbers.

Your 1/2^k point helps me understand why == does seem to work correctly
sometimes.

I wonder if we should be suspicious of >=. Imagine the horror if a= w/x >
b=y/z in fractions, but digitally a < b. Blech. Can that happen?


But, change the example's divisor from 40 to 30 [the fractional part
from 1/2 to 2/3]:

> (100*23)/30
[1] 76.671404

> 100*(23/30)
[1] 76.671404

Now the two operations give the same answer to the full precision
available.  So, it isn't "generally true true in R that (100*x)/y is
more accurate than 100*(x/y), if x > y."


The good news here is that round() gives same answer in both cases:)

I am looking for a case where the first method is less accurate than the
second. I expect that multiplying integers before dividing is never less
accurate. Sometimes it is more accurate.
`
Following your 1/2^k insight, you see why multiplying first is helpful in
some cases. Question is will situation get worse.

But Bert is right. I have to read more books.

I studied Golub and van Loan and came away with healthy fear of matrix
inversion. But when you look at user contributed regression packages, what
do you find? Matrix inversion and lots of X'X.


Paul Johnson
University of Kansask


The key (in your example) is a property of the way that floating point
arithmetic is implemented.


---JRG



On 04/21/2017 08:19 AM, Paul Johnson wrote:
> We all agree it is a problem with digital computing, not unique to R. I
> don't think that is the right place to stop.
>
> What to do? The round example arose in a real funded project where 2 R
> programs differed in results and cause was  that one person got 57 and
> another got 58. The explanation was found, but its less clear how to
> prevent similar in future. Guidelines, anyone?
>
> So far, these are my guidelines.
>
> 1. Insert L on numbers to signal that you really mean INTEGER. In R,
> forgetting the L in a single number will usually promote whole calculation
> to floats.
> 2. S3 variables are called 'numeric' if they are integer or double
storage.
> So avoid "is.numeric" and prefer "is.double".
> 3. == is a total fail on floats
> 4. Run print with digits=20 so we can see the less rounded number. Perhaps
> start sessions with "options(digits=20)"
> 5. all.equal does what it promises, but one must be cautious.
>
> Are there math habits we should follow?
>
> For example, Is it generally true in R that (100*x)/y is more accurate
than
> 100*(x/y), if x > y?   (If that is generally true, couldn't the R
> interpreter do it for the user?)
>
> I've seen this problem before. In later editions of the game theory
program
> Gambit, extraordinary effort was taken to keep values symbolically as
> integers as long as possible. Avoid division until the last steps. Same in
> Swarm simulations. Gary Polhill wrote an essay about the Ghost in the
> Machine along those lines, showing accidents from trusting floats.
>
> I wonder now if all uses of > or < with numeric variables are suspect.
>
> Oh well. If everybody posts their advice, I will write a summary.
>
> Paul Johnson
> University of Kansas
>
> On Apr 21, 2017 12:02 AM, "PIKAL Petr" <petr.pi...@precheza.cz> wrote:
>
>> Hi
>>
>> The problem is that people using Excel or probably other such
spreadsheets
>> do not encounter this behaviour as Excel silently rounds all your
>> calculations and makes approximate comparison without telling it does so.
>> Therefore most people usually do not have any knowledge of floating point
>> numbers representation.
>>
>>  Cheers
>> Petr
>>


_

[R] workaround for "package building problem Centos 7 & R 3.4.1"

[Paul Johnson]

This is for the problem I posted about last Friday.

First, the happy part, a workaround:

$ cd ~/R/x86_64-redhat-linux-gnu-library/3.4
$ ln -sf /usr/share/R/library/* .

After that, all of the packages are found by R CMD check.  R CMD check
looks in the ~/R/x86_64-redhat-linux-gnu-library/3.4 folder for
packages, but not in "/usr/share/R/library" (for whatever reason, I do
not know yet).

Why packages are not found in /usr/share/R/library is the unsolved
mystery. I don't see any difference in configuration between Centos
and Debian. In  Centos /usr/lib64/R/etc/Renviron we have the EPEL
version of the settings that Dirk Eddelbuettel worked out with R Core
for Debian packaging way back in the beginning of, well, the century:

R_LIBS_USER=${R_LIBS_USER-'~/R/x86_64-redhat-linux-gnu-library/3.4'}
R_LIBS_SITE=${R_LIBS_SITE-'/usr/local/lib/R/site-library:/usr/local/lib/R/library:/usr/lib64/R/library:/usr/share/R/library'}

R_LIBS_SITE looks good in that file. I know (for sure) these settings
were put in as adaptation to Dirk's Debian package because I sent them
to Spot Calloway back in the early days of the EPEL package.

I've found a similar report about this in the EPEL R packaging bugzilla:

https://bugzilla.redhat.com/show_bug.cgi?id=1457395

Appears there that the RPM package builds, even though the package
does not survive the check. Interesting:)

Summary of the problem:

R CMD check has unexpected .libPaths().  It does not use the
R_LIBS_SITE settings in Renviron.

The error message from R CMD check without "--as-cran" is informative,
as it plainly states that .libPaths() does not seem to include the
location where the packages are installed:

###
$ R CMD check kutils_1.19.tar.gz
* using log directory ‘/home/pauljohn/GIT/CRMDA/kutils/package/kutils.Rcheck’
* using R version 3.4.1 (2017-06-30)
* using platform: x86_64-redhat-linux-gnu (64-bit)
* using session charset: UTF-8
* checking for file ‘kutils/DESCRIPTION’ ... OK
* checking extension type ... Package
* this is package ‘kutils’ version ‘1.19’
* checking package namespace information ... OK
* checking package dependencies ... OK
* checking if this is a source package ... OK
* checking if there is a namespace ... OK
* checking for executable files ... OK
* checking for hidden files and directories ... OK
* checking for portable file names ... OK
* checking for sufficient/correct file permissions ... OK
* checking whether package ‘kutils’ can be installed ... OK
* checking installed package size ... OK
* checking package directory ... OK
* checking ‘build’ directory ... OK
* checking DESCRIPTION meta-information ... OK
* checking top-level files ... OK
* checking for left-over files ... OK
* checking index information ... OK
* checking package subdirectories ... OK
* checking R files for non-ASCII characters ... OK
* checking R files for syntax errors ... OK
* checking whether the package can be loaded ... OK
* checking whether the package can be loaded with stated dependencies ... OK
* checking whether the package can be unloaded cleanly ... OK
* checking whether the namespace can be loaded with stated dependencies ... OK
* checking whether the namespace can be unloaded cleanly ... OK
* checking loading without being on the library search path ... WARNING
Error: package or namespace load failed for ‘kutils’ in
loadNamespace(j <- i[[1L]], c(lib.loc, .libPaths()), versionCheck =
vI[[j]]):
 there is no package called ‘openxlsx’
Execution halted

It looks like this package has a loading problem when not on .libPaths:
see the messages for details.
* checking dependencies in R code ... OK
##

However, I promise all of the packages are installed, in "/usr/share/R/library".

Interactive R sessions do find the required packages, and they also
seem to be found by R programs run with R CMD BATCH and R with or
without "--vanilla".

$ R --vanilla -e ".libPaths()"

R version 3.4.1 (2017-06-30) -- "Single Candle"
Copyright (C) 2017 The R Foundation for Statistical Computing
Platform: x86_64-redhat-linux-gnu (64-bit)

> .libPaths()
[1] "/home/pauljohn/R/x86_64-redhat-linux-gnu-library/3.4"
[2] "/usr/lib64/R/library"
[3] "/usr/share/R/library"
>

I don't understand why the Renviron settings are found, actually, but
I'm counting my good fortune that they are.

Before I found the simple easy workaround, I went though the tedious
process of building and installing all of the packages in
~/R/x86_64-redhat-linux-gnu-library/3.4. That solves the problem as
well.  It seems possible to me that many people have packages
installed in that location anyway, so they are never bothered by the
problem in the first place.

Anyway, we have a workaround, but I don't understand what's wrong to
begin with or how to fix it the right way.

If you keep Googling, you see that this problem was experienced as
early as 2012.  It has happened to some of my personal R heroes,
actually. You know who you are :)

pj
-- 
Paul E. 

[R] package building problem Centos 7 & R 3.4.1

[Paul Johnson]

Dear everybody:

Packages that DO pass the package check on my Ubuntu 17.04 laptop with
R 3.4.1 and Mac OSX do not build on Centos 7 with R 3.4.1. I'm pretty
sure I have an environment defect, but cannot find it.

I find posts from various people about this problem since 2012.  But
I've checked the likely suspects. Among all the runs I pasted in
below, the most informative error message I found is the very last
one, which ends like this:

It looks like this package has a loading problem when not on .libPaths:
see the messages for details.

The messages don't really help me understand (all pasted in full
detail below), but I hope you'll say "aha", there's that environment
thing again.

I notice that package check proceeds differently without "--as-cran".
In some cases, package check passes. Other times, it fails
differently.

Here are some details, I hope you have guesses about where to check for trouble.

I suspected I had an error in my packages, so I picked some at random
from the CRAN list, things which depend on openxlsx. It turns out they
don't survive check either. Hopefully this makes it easy for you to
replicate (copy/paste :) )

$ wget http://rweb.quant.ku.edu/cran/src/contrib/reproducer_0.1.9.tar.gz
$ R CMD check --as-cran reproducer_0.1.9.tar.gz
* using log directory ‘/tmp/pj/reproducer.Rcheck’
* using R version 3.4.1 (2017-06-30)
* using platform: x86_64-redhat-linux-gnu (64-bit)
* using session charset: UTF-8
* using option ‘--as-cran’
* checking for file ‘reproducer/DESCRIPTION’ ... OK
* this is package ‘reproducer’ version ‘0.1.9’
* checking CRAN incoming feasibility ... WARNING
Maintainer: ‘Lech Madeyski ’

Insufficient package version (submitted: 0.1.9, existing: 0.1.9)
* checking package namespace information ... OK
* checking package dependencies ... OK
* checking if this is a source package ... OK
* checking if there is a namespace ... OK
* checking for executable files ... OK
* checking for hidden files and directories ... OK
* checking for portable file names ... OK
* checking for sufficient/correct file permissions ... OK
* checking whether package ‘reproducer’ can be installed ... OK
* checking installed package size ... OK
* checking package directory ... OK
* checking DESCRIPTION meta-information ... OK
* checking top-level files ... OK
* checking for left-over files ... OK
* checking index information ... OK
* checking package subdirectories ... OK
* checking R files for non-ASCII characters ... OK
* checking R files for syntax errors ... OK
* checking whether the package can be loaded ... OK
* checking whether the package can be loaded with stated dependencies ... OK
* checking whether the package can be unloaded cleanly ... OK
* checking whether the namespace can be loaded with stated dependencies ... OK
* checking whether the namespace can be unloaded cleanly ... OK
* checking loading without being on the library search path ... OK
* checking use of S3 registration ... OK
* checking dependencies in R code ... OK
* checking S3 generic/method consistency ... OK
* checking replacement functions ... OK
* checking foreign function calls ... OK
* checking R code for possible problems ... OK
* checking Rd files ... OK
* checking Rd metadata ... OK
* checking Rd line widths ... OK
* checking Rd cross-references ... OK
* checking for missing documentation entries ... OK
* checking for code/documentation mismatches ... OK
* checking Rd \usage sections ... OK
* checking Rd contents ... OK
* checking for unstated dependencies in examples ... OK
* checking contents of ‘data’ directory ... OK
* checking data for non-ASCII characters ... OK
* checking data for ASCII and uncompressed saves ... OK
* checking examples ... ERROR
Running examples in ‘reproducer-Ex.R’ failed
The error most likely occurred in:

> base::assign(".ptime", proc.time(), pos = "CheckExEnv")
> ### Name: boxplotAndDensityCurveOnHistogram
> ### Title: boxplotAndDensityCurveOnHistogram
> ### Aliases: boxplotAndDensityCurveOnHistogram
>
> ### ** Examples
>
> library(ggplot2)
Error in library(ggplot2) : there is no package called ‘ggplot2’
Execution halted
* checking for unstated dependencies in ‘tests’ ... OK
* checking tests ...
  Running ‘testthat.R’
 ERROR
Running the tests in ‘tests/testthat.R’ failed.
Complete output:
  > library(testthat)
  Error: package or namespace load failed for 'testthat' in
loadNamespace(j <- i[[1L]], c(lib.loc, .libPaths()), versionCheck =
vI[[j]]):
   there is no package called 'R6'
  Execution halted
* checking PDF version of manual ... OK
* DONE

Status: 2 ERRORs, 1 WARNING
See
  ‘/tmp/pj/reproducer.Rcheck/00check.log’
for details.

Contents of 00check.log are same:

The error most likely occurred in:
> base::assign(".ptime", proc.time(), pos = "CheckExEnv")
> ### Name: boxplotAndDensityCurveOnHistogram
> ### Title: boxplotAndDensityCurveOnHistogram
> ### Aliases: boxplotAndDensityCurveOnHistogram
>
> ### ** Examples
>
> library(ggplot2)
Error in library(ggplot2) : 

Re: [R] package building problem Centos 7 & R 3.4.1

[Paul Johnson]

On Aug 12, 2017 11:58 AM, "José Abílio Matos" <jaoma...@gmail.com> wrote:

On Friday, 11 August 2017 22.51.12 WEST Paul Johnson wrote:
> Dear everybody:
>
> Packages that DO pass the package check on my Ubuntu 17.04 laptop with
> R 3.4.1 and Mac OSX do not build on Centos 7 with R 3.4.1. I'm pretty
> sure I have an environment defect, but cannot find it.

Hi Paul,
the first error that you get is that ggplot2 is not installed:
Error in library(ggplot2) : there is no package called ‘ggplot2’

IIRC that package is not available in the EPEL repositories.


That's not it.
All of the concerned packages are installed and run fine. I demonstrated
that in original post. EPEL has nothing to do with this, install.packages
works fine.

PJ


I think that most of the other errors are similar.

Regards,
--
José Matos

[[alternative HTML version deleted]]

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

[[alternative HTML version deleted]]

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

Re: [R] Plotting log transformed predicted values from lme

[Paul Johnson]

In the rockchalk package, I have a function called newdata that will help
with this. Plenty of examples. Probably my predictOmatic function will just
work. Motivation is in the vignette.

Paul Johnson
University of Kansas

On Aug 9, 2017 11:23 AM, "Alina Vodonos Zilberg" <alina.vodo...@gmail.com>
wrote:

> Hi,
>
> I am performing meta-regression using linear mixed-effect model with the
> lme() function  that has two fixed effect variables;one as a log
> transformed variable (x)  and one as factor (y) variable, and two nested
> random intercept terms.
>
> I want to save the predicted values from that model and show the log curve
> in a plot ; predicted~log(x)
>
> mod<-lme(B~log(x)+as.factor(y), random=~1|cohort/Study,
> weights=varFixed(~I(SE^2)), na.action=na.omit, data=subset(meta),
>   control = lmeControl(sigma = 1, apVar = FALSE))
> summary(mod)
>
> newdat <- data.frame(x=seq(min(meta$x), max(meta$x),,118))  # I have 118
> observations. #How do I add the factor variable to my newdat?
> newdat$pred <- predict(mod, newdat,level = 0,type="response")
>
> plot(B ~ x, data=meta)
> lines(B ~ x, data=newdat)
>
> Can you please assist me ?
>
> Thank you!
>
> Alina
>
> [[alternative HTML version deleted]]
>
> __
> R-help@r-project.org mailing list -- To UNSUBSCRIBE and more, see
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide http://www.R-project.org/
> posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
>

[[alternative HTML version deleted]]

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

Re: [R] Install package "diagram"

[Paul Johnson]

On Wednesday, August 16, 2017, AbouEl-Makarim Aboueissa <
abouelmakarim1...@gmail.com> wrote:

> Dear All:
>
> I am trying to install the package "diagram". It is in the list. But when I
> selected the package to install it, it says:
>
> Question: "would you like to use a personal library instead?"
>
> You should say yes here. Let go into directory where you are allowed to
write.

 Other one fails because you are not the administrator.
PJ

> I selected No.
>
> Then it says
>
> > utils:::menuInstallPkgs()
> Warning in install.packages(NULL, .libPaths()[1L], dependencies = NA, type
> = type) :
>   'lib = "C:/Program Files/R/R-3.4.1/library"' is not writable
> Error in install.packages(NULL, .libPaths()[1L], dependencies = NA, type =
> type) :
>   unable to install packages
> >
>
>
> Any help will be appreciated.
>
>
> with many thanks
> abou
> --
> __
> AbouEl-Makarim Aboueissa, PhD
> 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.
>


-- 
Paul E. Johnson   http://pj.freefaculty.org
Director, Center for Research Methods and Data Analysis http://crmda.ku.edu

To write to me directly, please address me at pauljohn at ku.edu.

[[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] function pointers?

[Paul Johnson]

We have a project that calls for the creation of a list of many
distribution objects.  Distributions can be of various types, with
various parameters, but we ran into some problems. I started testing
on a simple list of rnorm-based objects.

I was a little surprised at the RAM storage requirements, here's an example:

N <- 1
closureList <- vector("list", N)
nsize = sample(x = 1:100, size = N, replace = TRUE)
for (i in seq_along(nsize)){
closureList[[i]] <- list(func = rnorm, n = nsize[i])
}
format(object.size(closureList), units = "Mb")

Output says
22.4 MB

I noticed that if I do not name the objects in the list, then the
storage drops to 19.9 MB.

That seemed like a lot of storage for a function's name. Why so much?
My colleagues think the RAM use is high because this is a closure
(hence closureList).  I can't even convince myself it actually is a
closure. The R source has

rnorm <- function(n, mean=0, sd=1) .Call(C_rnorm, n, mean, sd)

The storage holding 1 copies of rnorm, but we really only need 1,
which we can use in the objects.

Thinking of this like C,  I am looking to pass in a pointer to the
function.  I found my way to the idea of putting a function in an
environment in order to pass it by reference:

rnormPointer <- function(inputValue1, inputValue2){
object <- new.env(parent=globalenv())
object$distr <- inputValue1
object$n <- inputValue2
class(object) <- 'pointer'
object
}

## Experiment with that
gg <- rnormPointer(rnorm, 33)
gg$distr(gg$n)

ptrList <- vector("list", N)
for(i in seq_along(nsize)) {
ptrList[[i]] <- rnormPointer(rnorm, nsize[i])
}
format(object.size(ptrList), units = "Mb")

The required storage is reduced to 2.6 Mb. Thats 1/10 of the RAM
required for closureList.  This thing works in the way I expect

## can pass in the unnamed arguments for n, mean and sd here
ptrList[[1]]$distr(33, 100, 10)
## Or the named arguments
ptrList[[1]]$distr(1, sd = 100)

This environment trick mostly works, so far as I can see, but I have
these questions.

1. Is the object.size() return accurate for ptrList?  Do I really
reduce storage to that amount, or is the required storage someplace
else (in the new environment) that is not included in object.size()?

2. Am I running with scissors here? Unexpected bad things await?

3. Why is the storage for closureList so great? It looks to me like
rnorm is just this little thing:

function (n, mean = 0, sd = 1)
.Call(C_rnorm, n, mean, sd)


4. Could I learn (you show me?) to store the bytecode address as a
thing and use it in the objects?  I'd guess that is the fastest
possible way. In an Objective-C problem in the olden days, we found
the method-lookup was a major slowdown and one of the programmers
showed us how to save the lookup and use it over and over.

pj



-- 
Paul E. Johnson   http://pj.freefaculty.org
Director, Center for Research Methods and Data Analysis http://crmda.ku.edu

To write to me directly, please address me at pauljohn at ku.edu.

__
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] [ESS] M-x R gives no choice of starting dir

[Paul Johnson]

If Emacs is not asking for starting directory, it is very likely your
init file has this somewhere:

(setq ess-ask-for-ess-directory nil)


On Mon, Sep 11, 2017 at 3:23 PM, Enrico Schumann  
wrote:
> On Mon, 11 Sep 2017, Christian writes:
>
>> Hi,
>>
>> I experienced a sudden change in the behavior of M-x R in not giving
>> me the choice where to start R. May be that I botched my
>> preferences. I am using Aquamacs 3.3 on MacOS 10.12.6
>>
>> Christian
>
> I suppose you are using ESS? There is a variable called
> 'ess-ask-for-ess-directory', which controls whether
> M-x R prompts for a directory. Perhaps you have set
> this to nil?
>
> I also Cc the ESS-help mailing list and suggest that
> follow-up be sent there.
>
>
> Kind regards
>  Enrico
>
> --
> Enrico Schumann
> Lucerne, Switzerland
> http://enricoschumann.net
>
> __
> ess-h...@r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/ess-help



-- 
Paul E. Johnson   http://pj.freefaculty.org
Director, Center for Research Methods and Data Analysis http://crmda.ku.edu

To write to me directly, please address me at pauljohn at ku.edu.

__
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] looking for formula parser that allows coefficients

[Paul Johnson]

Can you point me at any packages that allow users to write a
formula with coefficients?

I want to write a data simulator that has a matrix X with lots
of columns, and then users can generate predictive models
by entering a formula that uses some of the variables, allowing
interactions, like

y ~ 2 + 1.1 * x1 + 3 * x3 + 0.1 * x1:x3 + 0.2 * x2:x2

Currently, in the rockchalk package, I have a function simulates
data (genCorrelatedData2), but my interface to enter the beta
coefficients is poor.  I assumed user would always enter 0's as
place holder for the unused coefficients, and the intercept is
always first. The unnamed vector is too confusing.  I have them specify:

c(2, 1.1, 0, 3, 0, 0, 0.2, ...)

I the documentation I say (ridiculously) it is easy to figure out from
the examples, but it really isnt.
It function prints out the equation it thinks you intended, thats
minimum protection against user error, but still not very good:

dat <- genCorrelatedData2(N = 10, rho = 0.0,
  beta = c(1, 2, 1, 1, 0, 0.2, 0, 0, 0),
  means = c(0,0,0), sds = c(1,1,1), stde = 0)
[1] "The equation that was calculated was"
y = 1 + 2*x1 + 1*x2 + 1*x3
 + 0*x1*x1 + 0.2*x2*x1 + 0*x3*x1
 + 0*x1*x2 + 0*x2*x2 + 0*x3*x2
 + 0*x1*x3 + 0*x2*x3 + 0*x3*x3
 + N(0,0) random error

But still, it is not very good.

As I look at this now, I realize expect just the vech, not the whole vector
of all interaction terms, so it is even more difficult than I thought to get the
correct input.Hence, I'd like to let the user write a formula.

The alternative for the user interface is to have named coefficients.
I can more or less easily allow a named vector for beta

beta = c("(Intercept)" = 1, "x1" = 2, "x2" = 1, "x3" = 1, "x2:x1" = 0.1)

I could build a formula from that.  That's not too bad. But I still think
it would be cool to allow formula input.

Have you ever seen it done?
pj
-- 
Paul E. Johnson   http://pj.freefaculty.org
Director, Center for Research Methods and Data Analysis http://crmda.ku.edu

To write to me directly, please address me at pauljohn at ku.edu.

__
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] Transforming data for nice output table

[Paul Johnson]

On Mon, Aug 20, 2018 at 2:17 PM David Doyle  wrote:
>
> Hello everyone,
>
> I'm trying to generate tables of my data out of R for my report.
>
> My data is setup in the format as follows and the example can be found at:
> http://doylesdartden.com/R/ExampleData.csv
>
> LocationDateYear  GW_Elevation
> 127(I)5/14/2006 2006   752.46
> 119(I)5/14/2006 2006   774.67
> 127(I)6/11/2007 2007   752.06
> 119(I)6/11/2007 2007   775.57
>
> I would like to generate a table that showed
>
> LocationGW_Elevation 2006GW_Elevation 2007GW_Elevation xxx.
>
> 119(I)774.67  775.57
>   
> 127(I)752.46  752.06
>   
>   XX   XX
>
>  Any thoughts on how to transform the data so it would be in this format??
>
> Thank you for your time
>
> David Doyle

Dear David

I'd consider studying R's reshape function, it was intended exactly
for this purpose. No reason to adventure into any user-contributed
tidy places to get this done.

dta <- read.csv("http://doylesdartden.com/R/ExampleData.csv;)
dta <- dta[c("Location", "Year", "GW_Elevation")]
dta.wide <- reshape(dta, direction = "wide", idvar = "Location",
v.names = "GW_Elevation", timevar = "Year")
head(dta.wide)

  Location GW_Elevation.2006 GW_Elevation.2007 GW_Elevation.2008
1   127(I)752.46NA757.50
2   119(S)774.67778.76776.40
3   132(I)759.45761.68764.27
4   132(S)761.77761.04765.44
5   111(I)753.52763.24764.24
6   111(S)766.18772.84767.41
  GW_Elevation.2009 GW_Elevation.2010 GW_Elevation.2011 GW_Elevation.2012
1759.90756.40759.05759.31
2777.59777.45778.21778.13
3761.90764.03763.63763.99
4761.21763.12762.69759.57
5750.85764.37762.99763.90
6769.77767.88767.95767.19
  GW_Elevation.2013 GW_Elevation.2014 GW_Elevation.2015 GW_Elevation.2016
1756.07756.66757.72757.66
2778.88778.28775.16778.28
3761.22762.81762.36764.46
4763.19763.87761.94763.90
5764.42761.65764.02762.93
6770.20767.25767.74766.87

The main difference between this and your stated target is that your
target column names have spaces in them, which are forbidden in
column names of data frames. Here R used a period for joining strings.
You can override
that if you want to with the reshape function, but usually I'd let the periods
happen.

If you do want to replace period with spaces, it can be done, but you
break the warranty
on other uses of a data frame. (Could get rid of underscore after GW
in same way)

colnames(dta.wide) <- sub("Elevation.", "Elevation ",
colnames(dta.wide), fixed = TRUE)

I'd not try to use that wide frame for many other purposes because of
the spaces, but it works well if you want to make a pleasant table out
of it. For example, xtable is my favorite:

library(xtable)

xt <- xtable(dta.wide)
print(xt)

The latex from that prints out beautifully in a document. The print
method for xtable has a file parameter if you want to save the file.

Good Luck

pj



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



-- 
Paul E. Johnson   http://pj.freefaculty.org
Director, Center for Research Methods and Data Analysis http://crmda.ku.edu

To write to me directly, please address me at pauljohn at ku.edu.

__
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] lm equivalent of Welch-corrected t-test?

[Paul Johnson]

Long ago, when R's t.test had var.equal=TRUE by default, I wrote some
class notes showing that the result was equivalent to a one predictor
regression model.  Because t.test does not default to var.equal=TRUE
these days, I'm curious to know if there is a way to specify weights
in an lm to obtain the same result as the Welch-adjusted values
reported by t.test at the current time.  Is there a WLS equivalent
adjustment with lm?

Here's example code to show that lm is same as t.test with var.equal.
The signs come out differently, but otherwise the effect estimate,
standard error, t value are same:


set.seed(234234)
dat <- data.frame(x = gl(2, 50, labels = c("F", "M")))
dat$err <- rnorm(100, 0, 1)
dat$y <- ifelse(dat$x == "F", 40 + dat$err, 44 + dat$err)
m1 <- lm(y ~ x, dat)
summary(m1)
m1.t <- t.test(y ~ x, dat, var.equal = TRUE)
m1.t
## diff matches regression coefficient
(m1.t.effect <- diff(m1.t$estimate))
## standard error matches regression se
m1.t.stderr <- m1.t.effect/m1.t$statistic

If you run that, you see lm output:

Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept)  39.9456 0.1180  338.65   <2e-16 ***
xM3.9080 0.1668   23.43   <2e-16 ***
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Residual standard error: 0.8341 on 98 degrees of freedom
Multiple R-squared:  0.8485,Adjusted R-squared:  0.8469
F-statistic: 548.8 on 1 and 98 DF,  p-value: < 2.2e-16

and t.test:

> m1.t <- t.test(y ~ x, dat, var.equal = TRUE)
> m1.t

Two Sample t-test

data:  y by x
t = -23.427, df = 98, p-value < 2.2e-16
alternative hypothesis: true difference in means is not equal to 0
95 percent confidence interval:
 -4.239038 -3.576968
sample estimates:
mean in group F mean in group M
   39.9455843.85358

> (m1.t.effect <- diff(m1.t$estimate))
mean in group M
   3.908003
> m1.t.effect/m1.t$statistic
mean in group M
 -0.1668129




-- 
Paul E. Johnson   http://pj.freefaculty.org
Director, Center for Research Methods and Data Analysis http://crmda.ku.edu

To write to me directly, please address me at pauljohn at ku.edu.

__
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 to obtain standardized betas in multiply-imputed data

[Paul Johnson]

Greetings.
I would adjust approach  to calculate standardized estimates for each
imputed set. Then summarize them . The way you are doing it here implies
that standardization concept applies to model list, which seems doubtful.
The empirical std. dev. of the variables differs among imputed data sets,
after all.

I suppose I mean to say lm.beta is not intended to receive a list of
regressions. Put standardization in the with() work done on each imputed
set. I suspect it is as easy as putting lm.beta in there. If there is
trouble, I have a standardize function in the rockchalk package. Unlike
lm.beta, it actually standardizes variables and runs regression. lm.beta
resales coefficients instead.

Paul Johnson
University of Kansas

On Wed, Sep 26, 2018, 5:03 AM CHATTON Anne via R-help 
wrote:

> Dear all,
>
> I am having problems in obtaining standardized betas on a multiply-imputed
> data set. Here are the codes I used :
> imp = mice(data, 5, maxit=10, seed=42, print=FALSE)
> FitImp <- with(imp,lm(y ~ x1 + x2 + x3))
> Up to here everything is fine. But when I ask for the standardized
> coefficients of the multiply-imputed regressions using this command :
> sdBeta <- lm.beta(FitImp)
> I get the following error message:
> Error in b * sx : argument non numérique pour un opérateur binaire
>
> Can anyone help me with this please?
>
> Anne
>
> __
> R-help@r-project.org mailing list -- To UNSUBSCRIBE and more, see
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide
> http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
>

[[alternative HTML version deleted]]

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

Re: [R] ANOVA in R

[Paul Johnson]

We cannot read your message. Should post pure text, not html. Hm, my phone
now may post html, must try to stop. Your R code not legible. It seems to
be output? Lines all run together. I tried find articles you mention, but
"not found" resulted.

You should use aov() for fitting, then get post hoc comparisons.

In car package, Anova function will help. I may teach Anova soon, we'll see
if I have better answer then.

Paul Johnson
University of Kansas

On Wed, Oct 10, 2018, 1:14 AM Thanh Tran  wrote:

> Hi eveyone,
> I'm studying about variance (ANOVA) in R and have some questions to share.
> I read an article investigating the effect of factors (temperature, Asphalt
> content, Air voids, and sample thickness) on the hardness of asphalt
> concrete in the tensile test (abbreviated as Kic). Each condition was
> repeated four times (4 samples). In the paper, the authors used MINITAB to
> analyze Anova. The authors use "adjusted sums of squares" calculate the
> p-value I try to use ANOVA in R to analyze this data and get the result as
> shown in Figure 4. The results are different from the results in the
> article. Some papers say that in R, the default for ANOVA analysis is to
> use "sequential sums of squares" to calculate the p-value.
> So please help the following two questions: 1 / Introduction to code in R
> for anova analysis uses "adjusted sums of squares". The main part of the
> command in R / myself is as follows: > Tem = as.factor (temperature) > Ac =
> as.factor (AC) > Av = as.factor (AV) > Thick = as.factor (Thickness) >
> Twoway = lm (KIC ~ Tem + Ac + Av + Thick + Stamp + Ac + Stamp + Av + Stamp
> + Thick + Ac * Av + Ac * Thick + Av * Thick) > anova (twoway) 2/ When to
> use "sequential sums of squares" and when to use "adjusted sums of
> squares". Some papers recommend using the "oa.design
> <
> https://www.youtube.com/redirect?q=http%3A%2F%2Foa.design%2F_token=AaSAPDY-5UAsoHxN6BdwfyIJ7R98MTUzOTIxNDg2OUAxNTM5MTI4NDY5=comments
> >"
> function in R to check for "orthogonal" designs. If not, use "adjusted sums
> of squares". I am still vague about this command, so look forward to
> everyone's suggestion. If you could answer all two of my questions, I would
> be most grateful. Ps: I have added a CSV file and the paper for practicing
> R. http://www.mediafire.com/file/e5oe54p2c2wd4bc/Saha+research.csv
>
> http://www.mediafire.com/file/39jlf9h539y9mdz/Homothetic+behaviour+investigation+on+fracture+toughness+of+asphalt+mixtures+using+semicircular+bending+test.pdf
>
> [[alternative HTML version deleted]]
>
> __
> R-help@r-project.org mailing list -- To UNSUBSCRIBE and more, see
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide
> http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
>

[[alternative HTML version deleted]]

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

Re: [ESS] how to edit openbugs model file?

[Paul Johnson]

Thanks, everybody.  I'm taking the suggestion Martin offered, since
this thing is required to be named *.txt in the current OpenBUGS.

Solution

Put this at the end of the *.txt file:

# Local Variables:
# mode: R
# End:

That works.

This question was raised by an interesting problem. This might be a
useful a warning to others using Emacs and OpenBugs.  Punch line:
character encoding must be ASCII.

The student typed in a bugs model file on his system, using Emacs.
The default character encoding was set as UTF8. Then his model files
kept crashing OpenBugs and I could not understand why. It warned about
unfamiliar characters,  I kept checking to see that file encoding was
UTF8.

Doh! Then I remembered that OpenBUGS asks for ASCII files!  Symbols
like "*" and "~" are given different encoding in ASCII and UTF8. The
"~" looks exactly the same, whether ASCII or UTF8, but luckily for us,
that student's computer showed the UTF8 "*" as a highlighted character
that was in a different line position.

If you put the cursor on a character and run C-u C-x =, then you see a
nice display about the encoding. There's an especially helpful post
about it:
http://stackoverflow.com/questions/10500323/how-to-see-the-files-encoding-in-emacs.

>From within Emacs, I cannot replicate that problem that the student
generated. No matter how I set buffer-file-encoding-system, I get
ASCII whenever I type "*" or "~". I'm guessing that because the
student is from another country, and English is not his native
language, that his Locale is different somehow.  But I can't see
what's different about it. Interesting and frustrating!

pj

On Sat, Dec 3, 2016 at 8:40 AM, Sparapani, Rodney  wrote:
> Hi Paul:
>
> Use BUGS mode which is part of ESS…
> (require ‘ess-bugs-d)
>
> Then name your file model.bug
>
> Rodney
>
> __
> ESS-help@r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/ess-help



-- 
Paul E. Johnson   http://pj.freefaculty.org
Director, Center for Research Methods and Data Analysis http://crmda.ku.edu

To write to me directly, please address me at pauljohn at ku.edu.

__
ESS-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/ess-help

[ESS] working directory lost when R session closes

[Paul Johnson]

Do you notice this:

cd into a folder, say "~/tmp/project/R" and start emacs with a file in
there. The working directory correctly shows "~/tmp/project/R".

Then launch an R session. When you quit the R session, and start a new
R session, the working directory changed, it becomes "~/tmp/project".

This is just a little inconvenient if you get an R session with some
crap in it and you close it down to start fresh, but the WD is no
longer correct. Its necessary to close emacs and re-open the file.

I'll paste in the Emacs session transcript to show what I mean. The
only thing I do after the q() is hit the big blue R button:


> getwd()
[1] "/home/pauljohn/GIT/rockchalk/package/rockchalk/R"
> q()
Save workspace image? [y/n/c]: n

Process R finished at Fri Jul 14 09:30:04 2017


R version 3.4.1 (2017-06-30) -- "Single Candle"
Copyright (C) 2017 The R Foundation for Statistical Computing
Platform: x86_64-pc-linux-gnu (64-bit)

R is free software and comes with ABSOLUTELY NO WARRANTY.
You are welcome to redistribute it under certain conditions.
Type 'license()' or 'licence()' for distribution details.

  Natural language support but running in an English locale

R is a collaborative project with many contributors.
Type 'contributors()' for more information and
'citation()' on how to cite R or R packages in publications.

Type 'demo()' for some demos, 'help()' for on-line help, or
'help.start()' for an HTML browser interface to help.
Type 'q()' to quit R.

> > if(identical(getOption('pager'), file.path(R.home('bin'), 'pager'))) # 
> > rather take the ESS one
+   options(pager='cat')
> options(STERM='iESS', str.dendrogram.last="'", editor='emacsclient', 
> show.error.locations=TRUE)
> getwd()
[1] "/home/pauljohn/GIT/rockchalk/package/rockchalk"
>


-- 
Paul E. Johnson   http://pj.freefaculty.org
Director, Center for Research Methods and Data Analysis http://crmda.ku.edu

To write to me directly, please address me at pauljohn at ku.edu.

__
ESS-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/ess-help

Re: [ESS] Ever consider changing indentation on #?

[Paul Johnson]

On Thu, May 4, 2017 at 3:51 AM, Stephen Eglen <sj...@cam.ac.uk> wrote:
>
> hi Paul,
>
>> "what the Hell, the folding is borked!". This is more and more common
>> with the adoption of markdown-style commentary within R files.
>
> can you give an example of the markdown-style commentary?
>
>

Here is an R file that has "Rstudio sections" that begin with #.  I just tested
with their gui and the magic recipe for howto is written below.  It is
something like
Emacs org mode, I suppose. Code folding. This is an example of a
file that gets "destroyed" by Emacs/ESS indentation

## Paul Johnson
## 20170504

## Open this file in Rstudio, you'll see what I mean

# Here is a section from Code-> insert section 

x <- rnorm(100)
y <- rnorm(100)


## The effect of "Code -> Insert section" is to create a left aligned single #.
## Apparently, it follows with dots out to right edge.

# Here is a second section from Code -> insert section 

z <- function(a, b, c){
y <-3
y + 7
}

##  After this, in editor Rstudio, then the sections can be folded (hidden)


# Third section 

## How to create a section without using Rstudio menu?
## TO create this section, I did not use Rstudio menu. I insert
## 1. one #
## 2. a space
## 3. some words
## 4. Four dashes

## the # at begin and dashes at end is signal to R studio this is a section.


I understand I can fix my init.el go avoid this, but it would be nicer
for my purpose
if this were in ESS itself, because I keep telling students to use it
and they are put off/discouraged
by this.

I'm not an Rstudio user, but try to cooperate with some of them :)

pj

>> Until now, I've said "too bad" when they complain about their code,
>> but now I'm thinking we ought to consider changing the behavior of
>> Emacs/ESS.  Does it really really need to do # indentation that way?
>> Why?  I've always thought it is odd. And counterintuitive. In what
>> backwards world would one suppose # gets the most indent, while ###
>> gets none?. I don't see any benefit in that # indentation that way. I
>> do see a big benefit in the ## indent at the code level. I'm only
>> quibbling about #.
>
> just for background context, ESS follows the Emacs convention for
> comments (### for column 0, ## for same level, # for col 32 or more), that is
> used widely in elisp.  It is probably fair to say though that these
> conventions are not widely known outside of Emacs.
>
>> Would you ever consider changing the ESS default for R files so # is
>> flush left?  Could I offer you some cash to consider it?
>>
>> I understand I *could* learn code for my emacs.init to do that, but it
>> seems like a change that would be broadly beneficial to new Emacs &
>> ESS users.  I preach the message about Emacs to the students here, but
>> little quirks like this seem like an unnecessary hassle.
>
> Its likely to be fixable, but first it might be worth hearing Martin's
> views on this, as he may deal with differing commenting styles a lot in
> R-core.
>
> Stephen



-- 
Paul E. Johnson   http://pj.freefaculty.org
Director, Center for Research Methods and Data Analysis http://crmda.ku.edu

To write to me directly, please address me at pauljohn at ku.edu.

__
ESS-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/ess-help

[ESS] starting directory not detected correctly anymore

[Paul Johnson]

Greetings.

I want to have ESS know the working directory from shell pwd, I don't
want it to ask me.

I have a setting in my init.el that used to work that way:

;; start R in current working directory, don't let R ask user:
(setq ess-ask-for-ess-directory nil)

Recently, it has stopped working. The symptom of the problem is that
the R session working directory trims off the last element in the
path. This is in Ubuntu 17.04 with emacs 24.5.1 and ess 17.11.

When I open a file from command line, say in

$ cd ~/tmp/R/
$ emacs testme.R

and then I hit the bug blue icon to start R, then getwd() shows
"~/tmp". The last element in the path is lost.  Same happens if I
start R with "Alt-x R", so don't hate me for liking your pretty blue R
button.

I notice that M-x eshell gets it right, it opens a shell in the /R
directory (the correct one). Also Emacs File "Open Directory" also
gets it right. Its just the inferior ESS *R* session that doesn't get
it right.

If I remove that line from init.el, then the ESS process stops and
asks me what directory what I want and it always guesses correctly.

I never saw this happen before today when opening a pre-existing R
file from the shell.  I have seem similar in past if I have an Emacs
session open and close the R session and re-start a new R session
without closing Emacs.  That second instance almost always has lost
the "R" from the file path, and shows the parent directory.

Would somebody try it and tell me if I've just gone all the way off
the cliff toward crazy?

Suggestions welcome (except concerning cliffs), thanks as usual.

pj


-- 
Paul E. Johnson   http://pj.freefaculty.org
Director, Center for Research Methods and Data Analysis http://crmda.ku.edu

To write to me directly, please address me at pauljohn at ku.edu.

__
ESS-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/ess-help

Re: [ESS] starting directory not detected correctly anymore

[Paul Johnson]

I was running the dev version of ESS in Melpa. I cleared that out.

With Ubuntu Emacs 24.5, I still see the trouble. However, problem
disappeared with Emacs 25.1, which I notice is now available in the
standard Ubuntu repository. Rodney S wrote and said he did not see
problem in RedHat with Emacs 25.2, that's why I tried that.

I have tested this repeatedly with emacs25 and all is good.  With
emacs24, however, same problem exists.

For me, problem is solved by updating Emacs, but if you want some
debug info, point me at directions to get you what you need and I will
do that.

pj

On Sat, Dec 2, 2017 at 12:20 PM, Lionel Henry <lionel@gmail.com> wrote:
> Startup in the project/package root is experimental and isn't
> enabled by default.
>
> Is it possible you are using a dev version from a month or so ago?
> For instance you could have installed it with melpa. The released
> version or the current dev version should not exhibit that behaviour.
>
> Best,
> Lionel
>
>
>> On 2 déc. 2017, at 18:43, Vitalie Spinu <spinu...@gmail.com> wrote:
>>
>>
>> I think this is a consequence of the recent feature which sets the default
>> directory to the package directory. It looks like it treats your tmp 
>> directory
>> as a package. @Lionel?
>>
>>  Vitalie
>>
>>>> On Thu, Nov 30 2017 18:59, Paul Johnson wrote:
>>
>>> Greetings.
>>
>>> I want to have ESS know the working directory from shell pwd, I don't
>>> want it to ask me.
>>
>>> I have a setting in my init.el that used to work that way:
>>
>>> ;; start R in current working directory, don't let R ask user:
>>> (setq ess-ask-for-ess-directory nil)
>>
>>> Recently, it has stopped working. The symptom of the problem is that
>>> the R session working directory trims off the last element in the
>>> path. This is in Ubuntu 17.04 with emacs 24.5.1 and ess 17.11.
>>
>>> When I open a file from command line, say in
>>
>>> $ cd ~/tmp/R/
>>> $ emacs testme.R
>>
>>> and then I hit the bug blue icon to start R, then getwd() shows
>>> "~/tmp". The last element in the path is lost.  Same happens if I
>>> start R with "Alt-x R", so don't hate me for liking your pretty blue R
>>> button.
>>
>>> I notice that M-x eshell gets it right, it opens a shell in the /R
>>> directory (the correct one). Also Emacs File "Open Directory" also
>>> gets it right. Its just the inferior ESS *R* session that doesn't get
>>> it right.
>>
>>> If I remove that line from init.el, then the ESS process stops and
>>> asks me what directory what I want and it always guesses correctly.
>>
>>> I never saw this happen before today when opening a pre-existing R
>>> file from the shell.  I have seem similar in past if I have an Emacs
>>> session open and close the R session and re-start a new R session
>>> without closing Emacs.  That second instance almost always has lost
>>> the "R" from the file path, and shows the parent directory.
>>
>>> Would somebody try it and tell me if I've just gone all the way off
>>> the cliff toward crazy?
>>
>>> Suggestions welcome (except concerning cliffs), thanks as usual.
>>
>>> pj
>>
>> __
>> ESS-help@r-project.org mailing list
>> https://stat.ethz.ch/mailman/listinfo/ess-help
>



-- 
Paul E. Johnson   http://pj.freefaculty.org
Director, Center for Research Methods and Data Analysis http://crmda.ku.edu

To write to me directly, please address me at pauljohn at ku.edu.

__
ESS-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/ess-help

Re: [ESS] [R] M-x R gives no choice of starting dir

[Paul Johnson]

If Emacs is not asking for starting directory, it is very likely your
init file has this somewhere:

(setq ess-ask-for-ess-directory nil)


On Mon, Sep 11, 2017 at 3:23 PM, Enrico Schumann  
wrote:
> On Mon, 11 Sep 2017, Christian writes:
>
>> Hi,
>>
>> I experienced a sudden change in the behavior of M-x R in not giving
>> me the choice where to start R. May be that I botched my
>> preferences. I am using Aquamacs 3.3 on MacOS 10.12.6
>>
>> Christian
>
> I suppose you are using ESS? There is a variable called
> 'ess-ask-for-ess-directory', which controls whether
> M-x R prompts for a directory. Perhaps you have set
> this to nil?
>
> I also Cc the ESS-help mailing list and suggest that
> follow-up be sent there.
>
>
> Kind regards
>  Enrico
>
> --
> Enrico Schumann
> Lucerne, Switzerland
> http://enricoschumann.net
>
> __
> ESS-help@r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/ess-help



-- 
Paul E. Johnson   http://pj.freefaculty.org
Director, Center for Research Methods and Data Analysis http://crmda.ku.edu

To write to me directly, please address me at pauljohn at ku.edu.

__
ESS-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/ess-help

[ESS] polymode indents everything when editing Rmd file?

[Paul Johnson]

I wonder if I've got polymode set up correctly. I have Ubuntu 17.10
Emacs 25.2 Ess 17.11, it works well to edit R files. Today I installed
polymode from melpa,

When I open an Rmd file, I can visually see a difference.  The code
chunks are likely shaded areas. Also, some keywords within the R
chunks are highlighted.  "library" appears as light green, and the
quoted arguments within function calls are different colors.

However, while editing the file, it seems like editing features don't
work correctly.  In the chunks that are R code, indentation works as
expected. However, outside the chunks, polymode has no indentation
benefit. I expect things that should be flush left should stay there,
but tabbing causes indentation. For example, where a line starts with

# R data import

or

```{r}

A tab in either line, at any point within the line, causes 4 spaces before ```.

Lately, I've been trying to get Emacs to work with python and that has
been a massive pain.  Its possible I broke polymode by fiddingly with
that. So here's the init.el:

(when (>= emacs-major-version 24)
  (require 'package)
  (package-initialize)
  (add-to-list 'package-archives
   '("melpa" . "http://melpa.milkbox.net/packages/;) t)
  (add-to-list 'package-archives
   '("melpa-stable" . "https://stable.melpa.org/packages/;)))

(package-initialize)
(elpy-enable)


;; Org-mode with R doesn't work without this
;; http://orgmode.org/worg/org-contrib/babel/how-to-use-Org-Babel-for-R.html
(custom-set-variables
 ;; custom-set-variables was added by Custom.
 ;; If you edit it by hand, you could mess it up, so be careful.
 ;; Your init file should contain only one such instance.
 ;; If there is more than one, they won't work right.
 '(ansi-color-names-vector
   ["#2e3436" "#a4" "#4e9a06" "#c4a000" "#204a87" "#5c3566"
"#729fcf" "#ec"])
 '(column-number-mode t)
 '(cua-auto-tabify-rectangles nil)
 '(cua-mode t nil (cua-base))
 '(custom-enabled-themes (quote (whiteboard)))
 '(inferior-R-args "--no-restore-history --no-save")
 '(org-babel-load-languages (quote ((emacs-lisp . t) (R . t
 '(org-confirm-babel-evaluate nil)
 '(package-selected-packages (quote (polymode elpygen elpy)))
 '(show-paren-mode t))

(setq default-tab-width 4)

;; Section II. Keyboard and mouse customization
;; Mouse and cursor in the usual Mac/Windows way
(delete-selection-mode t)
;; 
http://stackoverflow.com/questions/13036155/how-to-to-combine-emacs-primary-clipboard-copy-and-paste-behavior-on-ms-windows
(setq select-active-regions t)
;; Trying to make mouse middle-click only paste from primary
;; X11 selection, not clipboard and kill ring:
(global-set-key [mouse-2] 'mouse-yank-primary)
;; highlight does not alter Kill ring:
(setq mouse-drag-copy-region nil)
;; windows style binding C-x, C-v, C-c, C-z:
(cua-mode t)
;; Don't tabify after rectangle commands:
(setq cua-auto-tabify-rectangles nil)

;; Section III. Programming conveniences:
;; light-up matching parens:
(show-paren-mode t)
;; turn on syntax highlighting:
(global-font-lock-mode t)
;; Auto fill is TRUE in text modes:
(setq text-mode-hook (quote (turn-on-auto-fill text-mode-hook-identify)))


;; Section IV. ESS Emacs Statistics

(require 'ess-site)
;; soft require: no error if package not found
(require 'ess-smart-underscore nil 'noerror)

(require 'poly-R)
(require 'poly-markdown)

(add-to-list 'auto-mode-alist '("\\.Rmd" . poly-markdown+r-mode))
;; start R in current working directory, don't let R ask user:
(setq ess-ask-for-ess-directory nil)

;; Change shortcut "run this line" to use Shift-Return
;; Suggested by Vitalie Spinu 2013-09-30 to co-exist with Windows Emacs
(eval-after-load "ess-mode"
 '(progn
   (define-key ess-mode-map [(control return)] nil)
   (define-key ess-mode-map [(shift return)] 'ess-eval-region-or-line-and-step))
)

;; Help in frames?
(setq ess-help-own-frame 'one)

;; I want the *R* process in its own frame
(setq inferior-ess-own-frame t)
(setq inferior-ess-same-window nil)

;; Problem in Windows about menus not displaying
;; https://stat.ethz.ch/pipermail/ess-help/2012-September/008207.html
;;(if (eq system-type 'windows-nt)
 (defun my-R-execute-options ()
(ess-command "options(menu.graphics = FALSE)\n"))
(add-hook 'ess-R-post-run-hook 'my-R-execute-options)
;;)

;; Put # flush left
;; See gmail thread in ESS group, this from Ista
(add-hook 'ess-mode-hook
  (lambda()
;; don't indent comments
(setq ess-indent-with-fancy-comments nil)
;; turn on outline mode
(setq-local outline-regexp "[# ]+")
(outline-minor-mode t)))

;; Section V. Frames oriented Emacs

;; Make files opened from the menu bar appear in their own
;; frames. Adapted from Emacs menu-bar.el
(defun menu-find-existing ()
  "Edit the existing file FILENAME."
  (interactive)
  (let* ((mustmatch (not (and (fboundp 'x-uses-old-gtk-dialog)
  (x-uses-old-gtk-dialog
 (filename (car 

Re: [ESS] polymode indents everything when editing Rmd file?

[Paul Johnson]

On Tue, Feb 6, 2018 at 11:31 AM, Ista Zahn <istaz...@gmail.com> wrote:
> Hi Paul,
>
> Why do you expect polymode to provide an indention benefit? "All" it
> does in the .Rmd case is give you a buffer with markdown-mode parts
> and R-mode parts. Indentation in the R-mode parts is handled by ESS,
> and indentation in markdown-mode parts is handled by markdown-mode.
>
> Or are you trying to say that indentation in the markdown portion of a
> polymode buffer is different than indentation in a plain markdown-mode
> buffer?
>

Ugh. I see.  It is not a polymode problem.  I just don't like
markdown-mode so much.

If anybody has tips for working with Rmd files, I'm all ears at this point.

Best
pj
> Best,
> Ista
>
> On Tue, Feb 6, 2018 at 11:24 AM, Paul Johnson <pauljoh...@gmail.com> wrote:
>> I wonder if I've got polymode set up correctly. I have Ubuntu 17.10
>> Emacs 25.2 Ess 17.11, it works well to edit R files. Today I installed
>> polymode from melpa,
>>
>> When I open an Rmd file, I can visually see a difference.  The code
>> chunks are likely shaded areas. Also, some keywords within the R
>> chunks are highlighted.  "library" appears as light green, and the
>> quoted arguments within function calls are different colors.
>>
>> However, while editing the file, it seems like editing features don't
>> work correctly.  In the chunks that are R code, indentation works as
>> expected. However, outside the chunks, polymode has no indentation
>> benefit. I expect things that should be flush left should stay there,
>> but tabbing causes indentation. For example, where a line starts with
>>
>> # R data import
>>
>> or
>>
>> ```{r}
>>
>> A tab in either line, at any point within the line, causes 4 spaces before 
>> ```.
>>
>> Lately, I've been trying to get Emacs to work with python and that has
>> been a massive pain.  Its possible I broke polymode by fiddingly with
>> that. So here's the init.el:
>>
>> (when (>= emacs-major-version 24)
>>   (require 'package)
>>   (package-initialize)
>>   (add-to-list 'package-archives
>>'("melpa" . "http://melpa.milkbox.net/packages/;) t)
>>   (add-to-list 'package-archives
>>'("melpa-stable" . "https://stable.melpa.org/packages/;)))
>>
>> (package-initialize)
>> (elpy-enable)
>>
>>
>> ;; Org-mode with R doesn't work without this
>> ;; http://orgmode.org/worg/org-contrib/babel/how-to-use-Org-Babel-for-R.html
>> (custom-set-variables
>>  ;; custom-set-variables was added by Custom.
>>  ;; If you edit it by hand, you could mess it up, so be careful.
>>  ;; Your init file should contain only one such instance.
>>  ;; If there is more than one, they won't work right.
>>  '(ansi-color-names-vector
>>["#2e3436" "#a4" "#4e9a06" "#c4a000" "#204a87" "#5c3566"
>> "#729fcf" "#ec"])
>>  '(column-number-mode t)
>>  '(cua-auto-tabify-rectangles nil)
>>  '(cua-mode t nil (cua-base))
>>  '(custom-enabled-themes (quote (whiteboard)))
>>  '(inferior-R-args "--no-restore-history --no-save")
>>  '(org-babel-load-languages (quote ((emacs-lisp . t) (R . t
>>  '(org-confirm-babel-evaluate nil)
>>  '(package-selected-packages (quote (polymode elpygen elpy)))
>>  '(show-paren-mode t))
>>
>> (setq default-tab-width 4)
>>
>> ;; Section II. Keyboard and mouse customization
>> ;; Mouse and cursor in the usual Mac/Windows way
>> (delete-selection-mode t)
>> ;; 
>> http://stackoverflow.com/questions/13036155/how-to-to-combine-emacs-primary-clipboard-copy-and-paste-behavior-on-ms-windows
>> (setq select-active-regions t)
>> ;; Trying to make mouse middle-click only paste from primary
>> ;; X11 selection, not clipboard and kill ring:
>> (global-set-key [mouse-2] 'mouse-yank-primary)
>> ;; highlight does not alter Kill ring:
>> (setq mouse-drag-copy-region nil)
>> ;; windows style binding C-x, C-v, C-c, C-z:
>> (cua-mode t)
>> ;; Don't tabify after rectangle commands:
>> (setq cua-auto-tabify-rectangles nil)
>>
>> ;; Section III. Programming conveniences:
>> ;; light-up matching parens:
>> (show-paren-mode t)
>> ;; turn on syntax highlighting:
>> (global-font-lock-mode t)
>> ;; Auto fill is TRUE in text modes:
>> (setq text-mode-hook (quote (turn-on-auto-fill text-mode-hook-identify)))
>>
>>
>> ;; Section IV. ESS Emacs Statistics
>>
>> (req

[ESS] re-mapping C-RET in 18.10, avoiding conflict with CUA mode

[Paul Johnson via ESS-help]

I'm running Ubuntu-18.10 and the new packaging for elpa-ess provides
ESS 18.10.  I still like CUA mode rectangles and so I want to remap
C-RET to S-RET.

5 years ago Vitalie S told me how to do that. The stanza in init.el was:

;; Change shortcut "run this line" to use Shift-Return
;; Suggested by Vitalie Spinu 2013-09-30 to co-exist with Windows Emacs
(eval-after-load "ess-mode"
 '(progn
   (define-key ess-mode-map [(control return)] nil)
   (define-key ess-mode-map [(shift return)] 'ess-eval-region-or-line-and-step))
)

In Ess 18.10-1 this no longer works, I thought I'd fix by revising:

(eval-after-load "ess-mode"
 '(progn
   (define-key ess-mode-map [(control return)] nil)
   (define-key ess-mode-map [(shift return)]
'ess-eval-region-or-line-visibly-and-step))
)

But it does not help.

Can you advise me?

-- 
Paul E. Johnson   http://pj.freefaculty.org
Director, Center for Research Methods and Data Analysis http://crmda.ku.edu

To write to me directly, please address me at pauljohn at ku.edu.

__
ESS-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/ess-help

Re: [ESS] re-mapping C-RET in 18.10, avoiding conflict with CUA mode

[Paul Johnson via ESS-help]

Thanks, everybody.  I confirm this works.

I am especially delighted to see my former undergraduate student Pavel
crop up to help the old guy get the job done!

PJ

On Thu, Oct 25, 2018 at 3:00 AM Vitalie Spinu  wrote:
>
>
> Hi Paul,
>
> Ess-mode.el changed to ess.el but we are still discussing the change. For now
> just replace "ess-mode" with "ess" and it should be fine. The change to
> *-visibility-* is needed.
>
> In the future release both ess-mode and ess will work.
>
>
>   Vitalie
>
> >> On Wed, Oct 24 2018 14:58, Paul Johnson via ESS-help wrote:
>
> > I'm running Ubuntu-18.10 and the new packaging for elpa-ess provides
> > ESS 18.10.  I still like CUA mode rectangles and so I want to remap
> > C-RET to S-RET.
>
> > 5 years ago Vitalie S told me how to do that. The stanza in init.el was:
>
> > ;; Change shortcut "run this line" to use Shift-Return
> > ;; Suggested by Vitalie Spinu 2013-09-30 to co-exist with Windows Emacs
> > (eval-after-load "ess-mode"
> >  '(progn
> >(define-key ess-mode-map [(control return)] nil)
> >(define-key ess-mode-map [(shift return)] 
> > 'ess-eval-region-or-line-and-step))
> > )
>
> > In Ess 18.10-1 this no longer works, I thought I'd fix by revising:
>
> > (eval-after-load "ess-mode"
> >  '(progn
> >(define-key ess-mode-map [(control return)] nil)
> >(define-key ess-mode-map [(shift return)]
> > 'ess-eval-region-or-line-visibly-and-step))
> > )
>
> > But it does not help.
>
> > Can you advise me?



-- 
Paul E. Johnson   http://pj.freefaculty.org
Director, Center for Research Methods and Data Analysis http://crmda.ku.edu

To write to me directly, please address me at pauljohn at ku.edu.

__
ESS-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/ess-help

[ESS] Run roxygen examples: how now?

[Paul Johnson via ESS-help]

In functions with roxygen examples, until very recently I could run the
examples even though they had the roxygen comment "##'" at the beginning.
Emacs/ESS would, literally, let me do "next" over and over again to run
example code.

Today I realize that does not happen anymore, Emacs/ESS is treating those
things like comments and the focus skips to the function definition.

Can I make it do the other way now?

I'm in Ubuntu, Emacs 25.2.2, ESS 18.10.2,

I may have broken a feature while trying to learn how to use python or Rmd
documents. My ~/.emacs/elpa directory has accumulated a lot of junk:
auto-complete-20170125.245
company-20190116.1133
dash-20180910.1856
deferred-20170901.1330
ein-20190117.359
elpy-20181228.1721
epl-20180205.2049
ess-smart-underscore-20190110.1750
find-file-in-project-20181217.246
flycheck-20190108.351
gnupg
highlight-indentation-20181204.839
htmlize-20180923.1829
ivy-20190116.1140
js2-mode-20180724.801
julia-mode-20190111.2044
markdown-mode-20181229.1430
material-theme-20171123.1840
math-symbol-lists-20190102.1831
pkg-info-20150517.1143
poly-markdown-20181010.2137
polymode-20190102.1910
poly-noweb-20190102.1938
poly-R-20190111.1908
popup-20160709.1429
py-autopep8-20160925.1052
pyvenv-20181228.1722
r-autoyas-20140101.1510
request-20181129.1138
s-20180406.808
simple-httpd-20190110.1505
skewer-mode-20180706.1807
websocket-20180423.16
yasnippet-20181013.2122
yasnippet-20181015.1212
yasnippet-snippets-20181211.2219

-- 
Paul E. Johnson   http://pj.freefaculty.org
Director, Center for Research Methods and Data Analysis http://crmda.ku.edu

To write to me directly, please address me at pauljohn at ku.edu.

[[alternative HTML version deleted]]

__
ESS-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/ess-help

Re: [ESS] Run roxygen examples: how now?

[Paul Johnson via ESS-help]

On Wed, Jan 30, 2019 at 4:37 PM Alex Branham  wrote:

>
> On Wed 30 Jan 2019 at 12:34, Paul Johnson via ESS-help <
> ess-help@r-project.org> wrote:
>
> > In functions with roxygen examples, until very recently I could run the
> > examples even though they had the roxygen comment "##'" at the beginning.
> > Emacs/ESS would, literally, let me do "next" over and over again to run
> > example code.
> >
> > Today I realize that does not happen anymore, Emacs/ESS is treating those
> > things like comments and the focus skips to the function definition.
> >
> > Can I make it do the other way now?
> >
> > I'm in Ubuntu, Emacs 25.2.2, ESS 18.10.2,
>
> Can you describe how you're trying to run them? C-RET seems to work fine
> for me.
>
> Thanks,
> Alex
>


Hi, Alex

In my init file I have S-RET as the replacement for C-RET, but I've tested
this both ways and get same behavior.  After removing the init file, when I
do C-RET, then the one line I'm on will run, but the focus jumps over the
rest of the example code lines and stops at the first line of code that is
not commented.  This happens whether an R file is edited in a package in
developer mode or if the file is in a separate directory.

If nobody else sees this, then I guess it has to be something peculiar
about Ubuntu's Emacs 25.2.2 ??

pj

-- 
Paul E. Johnson   http://pj.freefaculty.org
Director, Center for Research Methods and Data Analysis http://crmda.ku.edu

To write to me directly, please address me at pauljohn at ku.edu.

[[alternative HTML version deleted]]

__
ESS-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/ess-help

Re: [ESS] Run roxygen examples: how now?

[Paul Johnson via ESS-help]

On Mon, Feb 4, 2019 at 9:22 AM Alex Branham  wrote:

>
> > In my init file I have S-RET as the replacement for C-RET, but I've
> tested
> > this both ways and get same behavior.  After removing the init file,
> when I
> > do C-RET, then the one line I'm on will run, but the focus jumps over the
> > rest of the example code lines and stops at the first line of code that
> is
> > not commented.  This happens whether an R file is edited in a package in
> > developer mode or if the file is in a separate directory.
>
> Can you try with the latest ESS version (i.e. get it from MELPA rather
> than use 18.10)?[1]
>
> If that still doesn't work, you might have to do a bit of
> troubleshooting to find the issue since it's working for me. Here's a
> list of things to look at:
>
> - make sure ess-roxy-mode is on
> - make sure ess-roxy-remove-roxy-re is in ess-presend-filter-functions
> - make sure ess-roxy-re looks like it matches roxy comments
>
> > If nobody else sees this, then I guess it has to be something peculiar
> > about Ubuntu's Emacs 25.2.2 ??
>
> I guess that's possible, but it seems pretty unlikely.
>
> Footnotes:
> [1]  https://melpa.org/#/getting-started
>
> Interesting. Problem does not appear in ess-18.10.3snapshot

-- 
Paul E. Johnson   http://pj.freefaculty.org
Director, Center for Research Methods and Data Analysis http://crmda.ku.edu

To write to me directly, please address me at pauljohn at ku.edu.

[[alternative HTML version deleted]]

__
ESS-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/ess-help

[ESS] Roxygen Update is "off by one" in function

[Paul Johnson via ESS-help]

With Ubuntu Emacs 25.2 and ESS 18.10.3snapshot, In development mode within
a package, I just saw this weird thing happen.

I was in one function and used the menu ESS /Roxygen /Update Template and
then the parameter list was updated by insertion of the parameters before
this one in the file. That's a new thing, I've never seen it before.

Here's what I start with:

##' Draw standard error bar for discrete variables
##'
##' Used with plotSlopes if plotx is discrete
##'
##' @param x The estimate of the center point
##' @param lwr The lower confidence interval bound
##' @param upr The upper confidence interval bound
##' @param arrow.width Arrowhead length must be specified in inches. See
?arrows
##' @param width Thickness of shaded column
##' @param col Color for a bar
##' @param opacity 120 is default, that's partial see through.
##' @return NONE
##' @export
##' @author Paul Johnson
##'
se.bars <- function(x, y, lwr, upr, width=1, col = col.se, opacity = 120,
lwd = 1) {
h <- 0.35 * width #half width
q <- 0.15 * width #quarter width
d <- 0.05 * width #shaded area
iCol <- col2rgb(col)
bCol <- mapply(rgb, red = iCol[1,], green = iCol[2,],
   blue = iCol[3,], alpha = opacity, maxColorValue = 255)
### sCol: shade color
sCol <-  mapply(rgb, red = iCol[1,], green = iCol[2,],
blue = iCol[3,], alpha = opacity/3, maxColorValue = 255)
lCol <- mapply(rgb, red = iCol[1,], green = iCol[2,],
   blue = iCol[3,], alpha = min(255, 2*opacity),
maxColorValue = 255)
if (!is.null(lwr)){
polygon(c(x-d, x+d, x+d,  x-d), c(rep(lwr, 2), rep(upr, 2)), col =
sCol, border = bCol)
## PROBLEM: lwd widths are not linked to lwd parameter.
##  arrows(x0 = x, y0 = lwr, y1 = upr, code=3, angle=90,
length=arrow.width, col = lCol, cex=0.5, lwd = 2)
lines(c(x, x), c(lwr, upr), col = lCol, cex=0.5, lwd = 2, lend = 1,
ljoin = 1)
lines(c(x, x) + c(-q, q), c(lwr, lwr), col = lCol, cex = 0.5, lwd =
lwd, lend = 1, ljoin = 1, lmitre = 2)
lines(c(x, x) + c(-q, q), c(upr, upr), col = lCol, cex = 0.5, lwd =
lwd, lend = 1, ljoin = 1, lmitre = 2)
}
lines(c(x, x) + c(-h, h), c(y, y), col = lCol, lwd = lwd + 1)
NULL
}

After ESS / Roxygen /Update, I get all the arguments inserted that don't
belong with this function
##'
##' Used with plotSlopes if plotx is discrete
##'
##' @param newdf
##' @param olddf
##' @param plotx
##' @param modx
##' @param modxVals
##' @param interval
##' @param plotPoints
##' @param plotLegend
##' @param legendTitle
##' @param col Color for a bar
##' @param llwd
##' @param opacity 120 is default, that's partial see through.
##' @param ...
##' @param x The estimate of the center point
##' @param lwr The lower confidence interval bound
##' @param upr The upper confidence interval bound
##' @param arrow.width Arrowhead length must be specified in inches. See
?arrows
##' @param width Thickness of shaded column
##' @return NONE
##' @export
##' @author Paul Johnson
##'

Those newly inserted arguments come from the one before.

Between functions in the R file, I have NULL because that was recommended
for Roxygen authors in the old days.  Nothing else special.

Can you reproduce?

PJ

-- 
Paul E. Johnson   http://pj.freefaculty.org
Director, Center for Research Methods and Data Analysis http://crmda.ku.edu

To write to me directly, please address me at pauljohn at ku.edu.

[[alternative HTML version deleted]]

__
ESS-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/ess-help