Re: [R] multilevel binary and ordered regression models

2013-06-06 Thread Rune Haubo
On 6 June 2013 00:13, Xu Jun  wrote:
> Dear r-helpers,
>
> I have two questions on multilevel binary and ordered regression models,
> respectively:
>
> 1. Is there any r function (like lmer or glmer) to run multilevel ordered
> regression models?

Yes, package ordinal will fit such models.

Cheers,
Rune

__
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] lineplot.ci{sciplot}: Split x.factor label in two lines

2013-06-06 Thread Kumar Mainali
Some of the levels in my X factor has more than one word (e.g., North
America) which I would like to split into two lines so that my plot does
not get too wide. Any help is appreciated.

lineplot.CI(x.factor = Continent, response = PCC.PrsPts, group =
PointSource, data = master2, cex = 1.2,
 xlab = "", ylab = "Sensitivity", ylim = c(0.37,1.01), cex.lab = 0.9, x.leg
= 0.9, y.leg = 0.6, cex.leg=0.75, ncol=1, cex.axis=0.85,
 main = "",
col =  c (colorRampPalette(c("red", "orange"))  (2),
 colorRampPalette(c("blue", "green"))  (2),
colorRampPalette(c("purple", "black")) (2)),
 pch = c(19,19,18,1817,17))



-- 
Section of Integrative Biology
University of Texas at Austin
Austin, Texas 78712, USA

[[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] sample size determination - comparing means for different size samples

2013-06-06 Thread David Winsemius
This is a duplicate question, right?

On Jun 6, 2013, at 12:58 PM, Rebecca Greenblatt wrote:

> Looking to determine sample sizes for both my experimental and control
> groups (I want only a small portion of my participants in my experimental
> condition) in order to compare population means. I would be able to
> estimate standard deviation beforehand.
> 
> I'm using the bpower function from the Hmisc R package, specifically
> bsamsize to compare population proportions, but was wondering if there was
> an equivalent function for comparing means.
> 
> -- 
> Rebecca Greenblatt
> 
>   [[alternative HTML version deleted]]
> 

David Winsemius
Alameda, CA, USA

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


[R] error running mvabund package

2013-06-06 Thread andres.holz
Dear All,

This is my first post, and probably (and hence apologies that) my question
is very silly!
 
I'm having issues with a the mvabund package
(http://cran.r-project.org/web/packages/mvabund/index.html), 
and would be great to get some help!

Here is the code (and files are attached):

library(mvabund)
##visualizing data
florabund <- read.csv("CPL_floristics_abund_v1d.csv", header=TRUE)
vegtype <- read.csv("CPL_vegtype_all_v1b.csv", header=TRUE)
transect <- read.csv("CPL_transect_all_v1b.csv", header=TRUE)

cplflorabund <-mvabund(florabund, check.rows=FALSE, check.names=TRUE)
vegtype=as.factor(vegtype$vegtype)
plot(cplflorabund~vegtype) #there is a column called "treatment", and
another called "block"
plot.mvabund(cplflorabund)

##Fitting predictive models
transect=as.factor(transect$transect)
flormat=as.matrix(florabund)
abund.nb <-manyglm(flormat~transect*vegtype, family="nagative.binomial") 

predict(abund.nb, type= "response") 

##Checking Model Assumptions
plot(abund.nb) 

meanvar.plot(flormat~transect, col=as.numeric (vegtype)) #direct plot of
mean-variace relationship

anova(abund.nb, p.uni="adjusted") 

It's after running the last code line  [>anova(abund.nb, p.uni="adjusted")],
when I get thiis error:

Error in XvarIn[nterms - i, varseq > i + minterm] <- 0 : 
  (subscript) logical subscript too long

Any idea what I'm doing wrong?  I've tried shortening the length of all
files, but didn't matter. 

Thanks so much for your time!!

Cheers,
Andrés 

errorrunningmvabund.zip
  



--
View this message in context: 
http://r.789695.n4.nabble.com/error-running-mvabund-package-tp4668879.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] sample size determination - comparing means for different size samples

2013-06-06 Thread Rebecca Greenblatt
Looking to determine sample sizes for both my experimental and control
groups (I want only a small portion of my participants in my experimental
condition) in order to compare population means. I would be able to
estimate standard deviation beforehand.

I'm using the bpower function from the Hmisc R package, specifically
bsamsize to compare population proportions, but was wondering if there was
an equivalent function for comparing means.

-- 
Rebecca Greenblatt

[[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] Distance-based non-parametric ANOVA

2013-06-06 Thread Mikhail Umorin
Hello, 

I am comparing treatments by comparing within group to between group distances 
like 
described in 

MJ Anderson. 2001. A new method for non-parametric multivariate analysis of 
variance. 
Austral Ecology 26: 32 -- 46.

The idea is to find the ratio of within group sum of distance^2 to between 
group sum of 
distance^2. Then find the distribution of the statistic by permuting the sample.

Has anyone done this before in R code? Are there any packages (I did checked on 
CRAN - none). Code?

I want to make sure I won't be inventing a bicycle before I write my own 
package.

Comments? Suggestions? (I have never written an R package before)

Thank you for your time, 

Mikhail.

[[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] Not sure this is something R could do but it feels like it should be.

2013-06-06 Thread Jim Lemon
On 06/07/2013 01:03 AM, Polwart Calum (COUNTY DURHAM AND DARLINGTON NHS 
FOUNDATION TRUST) wrote:

Some colleagues nationally have developed a system which means they can pick 
the optimal sets of doses for a drug.  The system could apply to a number of 
drugs.  But the actual doses might vary.  To try and explain this in terms that 
the average Joe on the street might understand if you have some amoxicillin 
antibiotic for a chest infection the normal dose for an adult is 250 to 500mg 
increased to maybe 1000mg in severe cases.

For a child it is dosed from a liquid and people usually go from 62.5mg, 125mg 
to 250mg although you could measure any volume you wanted.

What this new method has developed is a means to pick the "right" standard 
doses so what above is 62.5, 125, 250, 500, 1000.  However the method they've used is 
really engineered about ensure the jump between doses is correct - you'll notice that the 
list above is a doubling up method.

But you can also have a doubling up method that went 50, 100, 200, 400, 800, 
1600  and pretty much as many as you can think of depending on your starting 
point and there is no scientific means to pick that starting point.  So 
colleagues have developed their rather more complex equivalent of the doubling 
method to determine the doses they need but they need to know if they should 
start at 40, 50, 62.5 or some other number.

Once they have the starting number they can calculate all the other doses.  I 
realise R can do that, and I realise using a loop of possible starting numbers 
it can build all those options.

Each patient then has a theoretical dose they should get lets say that's 
10mg/kg and you might treat patients from 5 to 120kg.  They are then looking to 
calculate the variance for each dose range so if we take the 50, 100, 200, 400 
model and said you'd give 50mg to anyone needing 0?? to 75mg 100mg to anyone 
needing 76 - 150mg etc... from there they are taking that range and saying 
that's a 31% overdose to a 33% underdose.  Then they want to find if there is a 
starting number which minimises the extent of under and overdosing...

Anyone know of an existing stats function in R that can easily do that and almost then 
report from some inputs a single number that is the "best fit"?

Calum


Hi Calum,
I can only answer from the perspective of someone who calculated doses 
of alcohol for experimental subjects many years ago. It was not possible 
to apply a linear function across the range due to a number of factors. 
One is that BAC, which was the target value, is dependent upon the 
proportion of the weight that represents the water compartment of the 
body. This varies with both weight (heavier people typically have a 
higher proportion of fat) and sex (women also tend to have slightly more 
fat). The real monkey wrench in the works was absorption rate, which 
often made nonsense of my calculations. This may not be as important in 
therapeutic drugs, for we were aiming at a specified BAC at a certain 
time after dosing rather than an average level. However, I suspect that 
many therapeutic drugs have a different dose by weight for children (we 
weren't dosing children) and choosing a starting point at the bottom of 
the range would almost certainly introduce a systematic error. My 
intuition would be to anchor the dosage rate in the middle of the scale 
and then extrapolate in both directions (adults only, of course).


Jim

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


Re: [R] SPlus script

2013-06-06 Thread Rolf Turner

On 07/06/13 03:19, Scott Raynaud wrote:

I actually had tried placing arguments in the call but it didn't work.   
However, I did
not think about writing it to a variable and printing.  That seems to have done 
the
trick.  Funny, I don't remember having to do that before, but that's not 
surprising.



If I remember correctly --- haven't used Splus for decades --- this is a 
difference

between Splus and R.

In R the output of a function is returned *invisibly* if that function 
is called

from within another function.  And source() is one such other function.

So if you have a script, say "melvin" with the single line:

sin(42)

and in R you execute

source("melvin")

you will see no output.  If in another script, say "clyde" you have the
single line

print(sin(42))

and in R you execute

source("clyde")

you will see

[1] -0.9165215

In Splus, IIRC, the print() call is unnecessary.  I.e. you would get the 
same

result by sourcing "melvin" and "clyde".

Current Splus users may correct me if I am wrong about this.

cheers,

Rolf Turner

__
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] generating a bar chart with two axis for co-linear variable

2013-06-06 Thread arun
HI,
Not sure if this is what you wanted.

pdf("BarplotsNew.pdf")
library(plotrix)
lst2<-lapply(seq_len(ncol(dat1)),function(i){
                    Ctdat<- table(dat1[,i])
                    Ctdat1<-(Ctdat/sum(Ctdat))*100
                    dat2<-data.frame(Ctdat,Ctdat1,stringsAsFactors=FALSE)[,-3]
                    colnames(dat2)[3]<-"Rel.Freq"
                    dat2[,1]<- as.numeric(as.character(dat2[,1]))
                    with(dat2,twoord.plot(Var1,Rel.Freq,Freq,
                        lylim=c(0,100),rylim=c(0,1000),
                        ylab="Relative Frequency",
                        rylab="Frequency",main=paste("Bar 
plot:",colnames(dat1)[i],sep=" "),
                        
type=c("bar","l"),lcol=2,rcol=4,xtickpos=Var1,xticklab=Var1))
                        })
 dev.off()

A.K.



- Original Message -
From: arun 
To: Sudha Krishnan 
Cc: R help 
Sent: Thursday, June 6, 2013 2:44 PM
Subject: Re: [R]  generating a bar chart with two axis for co-linear variable

HI,

May be this helps:
dat1<- read.table("sampledata.txt",header=TRUE,sep=",",stringsAsFactors=FALSE)
pdf("Barplots.pdf")
 lst1<-lapply(seq_len(ncol(dat1)),function(i) {Ctdat<- 
table(dat1[,i]);Ctdat1<-(Ctdat/sum(Ctdat))*100;barplot(Ctdat1,ylim=c(0,100),xlab=colnames(dat1)[i],ylab="Relative
 Frequency",main=paste("Barplot:",colnames(dat1)[i],sep=" "))})
dev.off()
A.K.




- Original Message -
From: Sudha Krishnan 
To: "r-help@r-project.org" 
Cc: 
Sent: Thursday, June 6, 2013 2:37 AM
Subject: [R]  generating a bar chart with two axis for co-linear variable



Hello Dimitris,



I was goggling for some help on Sensitivity vs 1-specificity and saw your link.



I hope you can be of help to me in one of the issue that I am facing in 
generating combo chart(bar chart and plot). I am a novice and have some 
difficulty in getting this logic correct.





I am give a dataset (I am attaching a sample dataset).



I am using a barplot() and passing values for percentage frequency and the 
corresponding variables. I am struck here, what my function does is only 
calculate the frequency for the listed variables and not the frequency 
percentage. Is there a method or a script with which I can pass the frequency 
percent and the related values as category columns for x axis?



I will attach the graphs that I have generated so that you can suggest the 
better way.



Sampledata - Sampledata.txt

What my function does to calculate the frequency with category names in X axis 
- 1.png

My requirement is to generate percentage frequency of the variable in y1 and 
not the frequency itself. 2.png (where x categories are missing)





Thanks,

Sudha Krishnan



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


BarplotsNew.pdf
Description: Adobe PDF document
__
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] highlighted a certain time period on multiple plots

2013-06-06 Thread R. Michael Weylandt
On Thu, Jun 6, 2013 at 6:43 AM, Ye Lin  wrote:
> Hey All,
>
> I have a dataset like this:
>
> DatedayVar1Var2Var3Obs1/1/2013Tue23411/2/2013Wed23521/3/2013Thu24631/4/2013
> Fri24741/5/2013Sat24.5851/6/2013Sun24.9961/7/2013Mon25.31071/8/2013Tue25.711
> 81/9/2013Wed26.11291/10/2013Thu26.513101/11/2013Fri26.914111/12/2013Sat27.3
> 15121/13/2013Sun27.716131/14/2013Mon28.117141/15/2013Tue28.518151/16/2013Wed
> 28.919161/17/2013Thu29.320171/18/2013Fri29.721181/19/2013Sat210.12219
> 1/20/2013Sun210.523201/21/2013Mon210.924211/22/2013Tue211.325221/23/2013Wed2
> 11.726231/24/2013Thu212.127241/25/2013Fri212.528251/26/2013Sat212.92926
> 1/27/2013Sun213.330271/28/2013Mon213.731281/29/2013Tue214.132291/30/2013Wed2
> 14.533301/31/2013Thu214.93431
> Here is the code I use to plot:
>
> par(mar=c(10, 0.5, 0.5, 0.5))
> par(mfrow=c(3,1))
> plot(Var1~Obs,data=dat,xaxt="n",xlab="")
> plot(Var2~Obs,data=dat,xaxt="n",xlab="")
> plot(Var3~Obs,data=dat,xaxt="n",xlab="")
> axis(1,at=dat$Obs,label=dat$Date)
> mtext(1,text="Date",line=2.5)
> axis(1,at=dat$Obs,label=dat$day,line=4)
> mtext(1,text="Day of week",line=7)
>
> How can I remove the extra white space between plots and emphasize time
> periods=weekend with dashed area?? I have attached original output and
> ideal output I am looking for.

Take a look at my "xtsExtra" package off of R-forge:

## Using your data as Jean arranged it.
library(xtsExtra)

dat2 <- xts(dat[,c(3,4,5)], as.Date(dat[,1], format = "%m/%d/%Y"))

plot(dat2, yax.loc = "left")

and look at

example(plot.xts)

for how to do shading / lines as desired.

Cheers,
MW

>
> Any suggestion to emphasize/highlight weekend periods is really appreciate!
>
> 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.
>

__
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] Revolutions Blog: May Roundup

2013-06-06 Thread David Smith
Revolution Analytics staff write about R every weekday at the Revolutions blog:
 http://blog.revolutionanalytics.com
and every month I post a summary of articles from the previous month
of particular interest to readers of r-help.

In case you missed them, here are some articles related to R from the
month of May:

Billions of geotagged Tweets create a beautiful map of the world when
plotted with the ggmap package: http://bit.ly/111ZQq4

A review of Ryan Sheftel's talk at R/Finance, on how he uses R on the
trading desk at Credit Suisse: http://bit.ly/111ZQq5 . Also, a quick
take on some other talks at the meeting: http://bit.ly/111ZQq6

Two new R books to be published in July: Applied Predictive Modeling
(Kuhn and Johnson), and Dynamic Documents with R and knitr (Xie)
http://bit.ly/111ZQq3

Some organizations using R in 2013: Google, The New York Times,
Facebook, FourSquare, Twitter, John Deere, Zillow, FDA, NOAA and CFPB:
http://bit.ly/111ZOOX

Sentiment analysis on the Enron emails using R and Infinit.e used to
reveal signs of employee anxiety: http://bit.ly/111ZOOW

The 7th R/Rmetrics workshop will take place in Switzerland June
30-July 4: http://bit.ly/111ZOOY

Highlights from the Milwaukee Workshop on R and Bioinformatics:
http://bit.ly/111ZOOZ

An R script to simulate and animate an escape from a Zombie horde:
http://bit.ly/111ZQq7

R 3.0.1 has been released: http://bit.ly/111ZOP1

The May edition of the Revolution Analytics newsletter, including R
training courses and a new gaming webinar on June 13:
http://bit.ly/111ZQq8

Social network analysis in R, presented at the New Frontiers in
Computing Conference: http://bit.ly/111ZOP0

My somewhat controversial comparison of the terms "Statistics", "Data
Science" and "Business Intelligence" (with video for context):
http://bit.ly/111ZOP3

An updated list of resources for R users, including the top 3
resources for R beginners: http://bit.ly/111ZOP2

Noam Ross's useful guide to speeding up R code: http://bit.ly/111ZQqc

Trevor Hastie presents on lasso and elastic-net regularization with
the glmnet package at the Orange County R User Group:
http://bit.ly/111ZQqb

Video replay of the "What's new in Revolution R Enterprise 6.2"
webinar: http://bit.ly/111ZOP4

Further critiques of a SAS white paper benchmarking SAS and R:
http://bit.ly/111ZQqd

Video replay of my Strata presentation, "Real-time Big Data Analytics:
>From Deployment to Production": http://bit.ly/111ZOP5

Naive Bayes modeling for big data with the RevoScaleR package:
http://bit.ly/111ZQqe

An analysis of CRAN package updates and the most prolific package
maintainers: http://bit.ly/111ZQqf

The New York Times again uses R to develop an on-line data
visualization, this time on NFL draft data: http://bit.ly/111ZOP6

Some non-R stories in the past month included: chaotic double
pendulums (http://bit.ly/111ZQqg), Game of Thrones family trees
(http://bit.ly/111ZOP7), top Big Data accounts on Twitter
(http://bit.ly/111ZOP8), StackExchange for Open Data
(http://bit.ly/111ZOP9), logical fallacies (http://bit.ly/111ZOPa), an
animation of every known meteorite (http://bit.ly/111ZQqh) and places
singles meet (http://bit.ly/111ZOPb).

Meeting times for local R user groups (http://bit.ly/eC5YQe) can be
found on the updated R Community Calendar at: http://bit.ly/bb3naW

If you're looking for more articles about R, you can find summaries
from previous months at http://blog.revolutionanalytics.com/roundups/.
Join the Revolution mailing list at
http://revolutionanalytics.com/newsletter to be alerted to new
articles on a monthly basis.

As always, thanks for the comments and please keep sending suggestions
to me at da...@revolutionanalytics.com . Don't forget you can also
follow the blog using an RSS reader, or by following me on Twitter
(I'm @revodavid).

Cheers,
# David

--
David M Smith 
VP of Marketing, Revolution Analytics  http://blog.revolutionanalytics.com
Tel: +1 (650) 646-9523 (Seattle WA, USA)
Twitter: @revodavid
We're hiring! www.revolutionanalytics.com/careers

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


Re: [R] how to fit a glm model with all possible interactions between explanatory variables

2013-06-06 Thread Bert Gunter
?formula

-- Bert

On Thu, Jun 6, 2013 at 12:52 PM, Jack Luo  wrote:
> Hi,
>
> I am trying to find a way to fit a glm model with all the possible
> interaction terms between different variables, without typing all the
> X1:X2, X1:X3,  in the formula, is there a way in R to do that?
>
> Thanks,
>
> -Jack
>
> [[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.



-- 

Bert Gunter
Genentech Nonclinical Biostatistics

Internal Contact Info:
Phone: 467-7374
Website:
http://pharmadevelopment.roche.com/index/pdb/pdb-functional-groups/pdb-biostatistics/pdb-ncb-home.htm

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


Re: [R] how to fit a glm model with all possible interactions between explanatory variables

2013-06-06 Thread David Winsemius

On Jun 6, 2013, at 12:52 PM, Jack Luo wrote:

> Hi,
> 
> I am trying to find a way to fit a glm model with all the possible
> interaction terms between different variables, without typing all the
> X1:X2, X1:X3,  in the formula, is there a way in R to do that?

The "*" operator does that:

  ~ X1 * X2 * X3

-- 
David Winsemius
Alameda, CA, USA

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


[R] how to fit a glm model with all possible interactions between explanatory variables

2013-06-06 Thread Jack Luo
Hi,

I am trying to find a way to fit a glm model with all the possible
interaction terms between different variables, without typing all the
X1:X2, X1:X3,  in the formula, is there a way in R to do that?

Thanks,

-Jack

[[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] highlighted a certain time period on multiple plots

2013-06-06 Thread Adams, Jean
It is helpful if you provide your example data in a format that is easy for
R-help readers to reproduce, for example, using the dput() function.  For
example,

dat <- structure(list(Date = c("1/1/2013", "1/2/2013", "1/3/2013",
"1/4/2013",
 "1/5/2013", "1/6/2013", "1/7/2013", "1/8/2013", "1/9/2013", "1/10/2013",
"1/11/2013", "1/12/2013", "1/13/2013", "1/14/2013", "1/15/2013",
 "1/16/2013", "1/17/2013", "1/18/2013", "1/19/2013", "1/20/2013",
"1/21/2013", "1/22/2013", "1/23/2013", "1/24/2013", "1/25/2013",
 "1/26/2013", "1/27/2013", "1/28/2013", "1/29/2013", "1/30/2013",
"1/31/2013"), day = c("Tue", "Wed", "Thu", "Fri", "Sat", "Sun",
 "Mon", "Tue", "Wed", "Thu", "Fri", "Sat", "Sun", "Mon", "Tue",
"Wed", "Thu", "Fri", "Sat", "Sun", "Mon", "Tue", "Wed", "Thu",
 "Fri", "Sat", "Sun", "Mon", "Tue", "Wed", "Thu"), Var1 = c(2L,
2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L,
 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L), Var2 = c(3,
3, 4, 4, 4.5, 4.9, 5.3, 5.7, 6.1, 6.5, 6.9, 7.3, 7.7, 8.1, 8.5,
 8.9, 9.3, 9.7, 10.1, 10.5, 10.9, 11.3, 11.7, 12.1, 12.5, 12.9,
13.3, 13.7, 14.1, 14.5, 14.9), Var3 = 4:34, Obs = 1:31), .Names = c("Date",
 "day", "Var1", "Var2", "Var3", "Obs"), class = "data.frame", row.names =
c(NA,
-31L))

Reduce the mar argument settings in the par() function to reduce the space
between plots, and use the oma argument to increase the space at the bottom
of the page.  You could use a for() loop to add vertical lines separately
to each plot.  Not quite what you had asked for, but very easy to implement.

par(mar=c(0.5, 5, 0.5, 0.5), oma=c(8, 0, 0, 0), mfrow=c(3, 1))
for(i in 1:3) {
colm <- match(paste0("Var", i), names(dat))
 plot(dat[, colm] ~ Obs, data=dat, xaxt="n", xlab="", ylab=paste0("Var", i))
abline(v=dat$Obs[dat$day %in% "Sat"]-0.5, lty=2, lwd=3, col="gray")
 abline(v=dat$Obs[dat$day %in% "Sun"]+0.5, lty=2, lwd=3, col="gray")
}
axis(1,at=dat$Obs,label=dat$Date)
mtext(1,text="Date",line=2.5)
axis(1,at=dat$Obs,label=dat$day,line=4)
mtext(1,text="Day of week",line=7)

Jean


On Thu, Jun 6, 2013 at 12:43 AM, Ye Lin  wrote:

> Hey All,
>
> I have a dataset like this:
>
> DatedayVar1Var2Var3Obs1/1/2013Tue23411/2/2013Wed23521/3/2013Thu24631/4/2013
>
> Fri24741/5/2013Sat24.5851/6/2013Sun24.9961/7/2013Mon25.31071/8/2013Tue25.711
> 81/9/2013Wed26.11291/10/2013Thu26.513101/11/2013Fri26.914111/12/2013Sat27.3
>
> 15121/13/2013Sun27.716131/14/2013Mon28.117141/15/2013Tue28.518151/16/2013Wed
> 28.919161/17/2013Thu29.320171/18/2013Fri29.721181/19/2013Sat210.12219
>
> 1/20/2013Sun210.523201/21/2013Mon210.924211/22/2013Tue211.325221/23/2013Wed2
> 11.726231/24/2013Thu212.127241/25/2013Fri212.528251/26/2013Sat212.92926
>
> 1/27/2013Sun213.330271/28/2013Mon213.731281/29/2013Tue214.132291/30/2013Wed2
> 14.533301/31/2013Thu214.93431
> Here is the code I use to plot:
>
> par(mar=c(10, 0.5, 0.5, 0.5))
> par(mfrow=c(3,1))
> plot(Var1~Obs,data=dat,xaxt="n",xlab="")
> plot(Var2~Obs,data=dat,xaxt="n",xlab="")
> plot(Var3~Obs,data=dat,xaxt="n",xlab="")
> axis(1,at=dat$Obs,label=dat$Date)
> mtext(1,text="Date",line=2.5)
> axis(1,at=dat$Obs,label=dat$day,line=4)
> mtext(1,text="Day of week",line=7)
>
> How can I remove the extra white space between plots and emphasize time
> periods=weekend with dashed area?? I have attached original output and
> ideal output I am looking for.
>
> Any suggestion to emphasize/highlight weekend periods is really appreciate!
>
> 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.
>
>

[[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] generating a bar chart with two axis for co-linear variable

2013-06-06 Thread arun
HI,

May be this helps:
dat1<- read.table("sampledata.txt",header=TRUE,sep=",",stringsAsFactors=FALSE)
pdf("Barplots.pdf")
 lst1<-lapply(seq_len(ncol(dat1)),function(i) {Ctdat<- 
table(dat1[,i]);Ctdat1<-(Ctdat/sum(Ctdat))*100;barplot(Ctdat1,ylim=c(0,100),xlab=colnames(dat1)[i],ylab="Relative
 Frequency",main=paste("Barplot:",colnames(dat1)[i],sep=" "))})
dev.off()
A.K.




- Original Message -
From: Sudha Krishnan 
To: "r-help@r-project.org" 
Cc: 
Sent: Thursday, June 6, 2013 2:37 AM
Subject: [R]  generating a bar chart with two axis for co-linear variable



Hello Dimitris,



I was goggling for some help on Sensitivity vs 1-specificity and saw your link.



I hope you can be of help to me in one of the issue that I am facing in 
generating combo chart(bar chart and plot). I am a novice and have some 
difficulty in getting this logic correct.





I am give a dataset (I am attaching a sample dataset).



I am using a barplot() and passing values for percentage frequency and the 
corresponding variables. I am struck here, what my function does is only 
calculate the frequency for the listed variables and not the frequency 
percentage. Is there a method or a script with which I can pass the frequency 
percent and the related values as category columns for x axis?



I will attach the graphs that I have generated so that you can suggest the 
better way.



Sampledata - Sampledata.txt

What my function does to calculate the frequency with category names in X axis 
- 1.png

My requirement is to generate percentage frequency of the variable in y1 and 
not the frequency itself. 2.png (where x categories are missing)





Thanks,

Sudha Krishnan



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


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


Re: [R] reshape2 issue

2013-06-06 Thread Adams, Jean
library(reshape2)
?cast

"Use acast or dcast depending on whether you want vector/matrix/array
output or data frame output."

Jean


On Wed, Jun 5, 2013 at 9:35 AM, Bruce Miller  wrote:

> Hi all,
>
> I am revisiting using reshape2 to aggregate critter (bats) occurrences
> by time blocks and have the final output for each row sum a percentage
> such that the row sums total 1 (100%).
>
> Reshape2 library is loaded on start.
>
> The old code I had stored in NotePad++ seems to no longer work as I have
> the following error:
>
>  > input.cast <- cast(input.melt, Species ~ Date, fun.aggregate = sum)
> *Error: could not find function "cast"*
>
> The code saved in NotePad++ is:
>
> require(reshape2)
> input.melt <- melt(input)
> input.cast <- cast(input.melt, Species ~ Date, fun.aggregate = sum)
> input.percent <- input.cast[, -1] / rowSums(input.cast[, -1]) * 100
> rownames(input.percent) <- input.cast[,1]
>
> The data is read in OK:
>  >Escameca_temporal <- read.table("C:/=Bat data
> working/=Nica_new/Escameca_temporal.CSV",header=T,sep=",",quote="")
> #then renamed : >input <-Escameca_temporal
> The variable view shows Species (6 letter text code) as Factor;Location
> code as integer; Date as Factor; Time as Factor.
>
> Data looks like this (sub sample only)
>
> Species LocationDateTime
> Cynmex  49375/3/201223:42
> Cynmex  49375/3/201223:43
> Cynmex  49375/3/201223:43
> Eptfur  49375/2/20120:10
> Eptfur  49375/2/20120:12
> Eptfur  49375/2/20120:15
> Eptfur  49375/2/20120:16
> Eptfur  49375/2/20120:21
> Eptfur  49375/2/20120:27
> Eptfur  49375/2/20120:27
> Eptfur  49375/2/20120:36
> Eptfur  49375/2/20120:37
> Eptfur  49375/2/20120:38
> Eptfur  49375/2/20120:42
> Eptfur  49375/2/20120:57
> Eptfur  49375/2/20120:58
> Eptfur  49375/2/20120:59
> Eptfur  49375/2/20121:01
> Eptfur  49375/2/20121:03
> Eptfur  49375/2/20121:08
> Eptfur  49375/4/20122:02
> Eptfur  49375/4/20122:03
> Eptfur  49375/4/20122:03
> Eptfur  49375/4/20122:04
> Eptfur  49375/4/20122:04
> Eptfur  49375/4/20122:08
> Eptfur  49241/8/20122:09
> Eptfur  49375/4/20122:09
> Eptfur  49241/8/20122:10
> Eptfur  49241/8/20122:10
> Eptfur  49241/8/20122:10
> Eptfur  49241/8/20122:10
> Eptfur  49241/8/20122:11
> Eptfur  49241/8/20122:11
> Eptfur  49241/8/20122:11
> Eptfur  49241/8/20122:11
> Eptfur  49241/8/20122:11
> Eptfur  49241/8/20122:12
> Eptfur  49241/8/20122:12
> Eptfur  49241/8/20122:12
> Eptfur  49241/8/20122:12
> Eptfur  49241/8/20122:12
> Eptfur  49241/8/20122:13
> Eptfur  49241/8/20122:13
> Eptfur  49375/4/20122:13
> Eptfur  49241/8/20122:14
> Eptfur  49241/8/20122:14
> Eptfur  49241/8/20122:14
> Eptfur  49241/8/20122:14
> Eptfur  49375/4/20122:14
> Eptfur  49241/8/20122:15
> Eptfur  49241/8/20122:15
> Eptfur  49241/8/20122:15
> Eptfur  49241/8/20122:16
> Eptfur  49241/8/20122:16
> Eptfur  49241/8/20122:16
> Eptfur  49241/8/20122:16
> Eptfur  49241/8/20122:17
>
>
> Any suggestions-help appreciated.
>
> Bruce
>
> [[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.


[R] Add prediction interval forest plot (package "meta")

2013-06-06 Thread Bernd Weiss
Dear all,

I am struggling to add a prediction interval to a forest plot that was
created with forest.meta(), package "meta".

I checked the source of forest.meta() and realized that it is heavily
relying on grid. I am lacking any experience with grid graphics. So, I
am having difficulties to find out where the random effects estimate is
actually plotted in order to add prediction intervals.

Any help is much appreciated!

Bernd


R version 3.0.1 Patched (2013-05-28 r62825)
Platform: x86_64-w64-mingw32/x64 (64-bit)

locale:
[1] LC_COLLATE=German_Germany.1252  LC_CTYPE=German_Germany.1252
LC_MONETARY=German_Germany.1252
[4] LC_NUMERIC=CLC_TIME=German_Germany.1252

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

other attached packages:
[1] meta_2.3-0

loaded via a namespace (and not attached):
[1] tools_3.0.1

__
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] OT: FAQ 7.31

2013-06-06 Thread Sarah Goslee
Or, how to tell thinking like humans from thinking like computers:

www.smbc-comics.com/index.php?db=comics&id=2999


-- 
Sarah Goslee
http://www.functionaldiversity.org

__
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] RBloomberg: unable to use 'Corp' to fetch data

2013-06-06 Thread Chirag Maru
I am trying to get new Cusip (registered security) from old .  This can be done 
in excel using BDP command.

For example:
=BDP("000361AJ  Corp","ID_CUSIP"), will give 000361AK1.

Can this be done using RBloomberg?  Any help will be much appreciated.

Chirag Maru
Senior Quantitative Analyst
IRON Financial
847-715-3221 Direct | 847-715-3321 Fax | 847-715-3200 Main
630 Dundee Road, Suite 200, Northbrook, IL 60062
chirag.m...@ironfinancial.com | 
www.ironfinancial.com




The information transmitted is intended solely for the individual or entity to 
which it is addressed and may contain confidential and/or privileged material. 
Any review, retransmission, dissemination or other use of or taking action in 
reliance upon this information by persons or entities other than the intended 
recipient is prohibited. If you have received this email in error please 
contact the sender and delete the material from any computer.

Any information contained herein is neither an offer to sell nor a solicitation 
to buy any interest in any investment fund. An offer can only be made by the 
approved offering memorandum, which contains important information concerning 
risk factors and other material information and must be read carefully before 
any decision to invest is made. Securities and derivatives trading are 
speculative and involve a risk of substantial loss. Past results are not 
necessarily indicative of future performance.

All incoming and outgoing e-mails are archived and may be reviewed and/or 
produced at the request of regulators or in connection with civil litigation. 
IRON Holdings, LLC. accepts no liability for any errors or omissions arising as 
a result of transmission.

[[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] Conditional coloring of area between curves

2013-06-06 Thread Duncan Murdoch

On 06/06/2013 10:41 AM, Roland Pape wrote:

Dear list,

I have two time series of temperatures from different sites plotted in the same 
diagram and would like to color the area between the curves differently, 
dependent on whether site 1 is cooler than site 2 (colored let's say in blue) 
or warmer (red). While polygone() works fine to color the enclosed area in 
general, I'm struggling with this conditional coloring. Any help is greatly 
appreciated!


Suppose the series are named "site1" and "site2", and they have common 
times in a variable "times".  Then the following should do what you want:


smaller <- pmin(site1, site2)
plot(x, site1, ylim = range(c(site1, site2)), type="n")
polygon(c(x, rev(x)), c(smaller, rev(site1)), col="red")
polygon(c(x, rev(x)), c(smaller, rev(site2)), col="blue")

If the times for the two series are different it's a little harder; 
first you need to give them common times, and that will depend on how 
you decide to evaluate the values between observations. Probably linear 
interpolation (using approx()) is fine, but it's up to you.


Duncan Murdoch

__
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] Optim seems not to work properly in foreach

2013-06-06 Thread Simon Zehnder
Hi Ilai,

after you sent this message I tried your code as well and it worked. As a 
result, I reconsidered the code written by me and of course also found the 
error in my function simulateEKOP. So for other people assuming errors or ill 
behavior in base functions: Forget it: these functions are tested this often, 
that the probabilitty is almost 1, the error is in your own code. 

Thanks for the help

Simon


On Jun 3, 2013, at 10:53 PM, ilai  wrote:

> On Mon, Jun 3, 2013 at 11:37 AM, Simon Zehnder 
> 
> ... [Some not minimal, self contained, reproducible code]...
>
> Data simulation and thecreation of startpar works fine, but the parameters in 
> res$par are always the start parameters. If I run the same commands directly 
> on the shell I get in res$par the optimized parameters - only inside the 
> foreach loop optim seems not to work. What could that be?
> 
> Don't know, but but this makes me doubt it has anything to do with optim 
> being inside foreach:   
> 
> fr <- function(x) { 
>  x1 <- x[1] ; x2 <- x[2] 
>  100 * (x2 - x1 * x1)^2 + (1 - x1)^2
>  }
> grr <- function(x) {
>  x1 <- x[1] ; x2 <- x[2]
>  c(-400 * x1 * (x2 - x1 * x1) - 2 * (1 - x1) , 200 * (x2 - x1 * x1))
>  }
> library("doMC") 
> registerDoMC(2)
> RNGkind("L'Ecuyer")
> set.seed(54321)
> foreach(i = 1:2) %do% {
>   ret <- foreach(j = 1:2) %do%{
>strtpar <- c(-2,2)+rnorm(2) 
>optim(strtpar, fr, grr, method = 
> "L-BFGS-B",control=list(trace=TRUE))$par
>   } 
>   ret 
> } 
> 
> Also, wouldn't you want to register 4 cores by default if nesting 2 loops of 
> 2 ? (to comment on the wisdom of doing so in terms of overhead is beyond my 
> expertise)  
>  
> HTH 
> 
> 
> 
> Best
> 
> Simon
> 
> __
> 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] generate simple function with pre-defined constants

2013-06-06 Thread Joshua Wiley
On Thu, Jun 6, 2013 at 9:05 AM, William Dunlap  wrote:
>
> I said the force was 'required' in the sense that without it
> the function will fail to do what you want in some situations.

With the Force on your side, functions always do what you want.

>
> Bill Dunlap
> Spotfire, TIBCO Software
> wdunlap tibco.com

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


Re: [R] generate simple function with pre-defined constants

2013-06-06 Thread Mark Leeds
gotcha. thanks.


On Thu, Jun 6, 2013 at 12:05 PM, William Dunlap  wrote:

>  I said the force was 'required' in the sense that without it
>
> the function will fail to do what you want in some situations.
>
> It doesn't make sense to write a function that you know will
>
> fail sometimes when you know an easy way to make it work
>
> in all situations.
>
> ** **
>
> Bill Dunlap
>
> Spotfire, TIBCO Software
>
> wdunlap tibco.com
>
> ** **
>
> *From:* Mark Leeds [mailto:marklee...@gmail.com]
> *Sent:* Thursday, June 06, 2013 8:57 AM
> *To:* William Dunlap
> *Cc:* Liviu Andronic; r-help@r-project.org Help
> *Subject:* Re: [R] generate simple function with pre-defined constants
>
> ** **
>
> hi bill: I understand what you're doing but, atleast for this case, I
> checked and you don't need the force this one. it works without it. so, I
> think the force requirement applies only when you're building them up with
> the lapply. but definitely I'm opened to clarification. thanks.
>
>
>
> 
>
> ** **
>
> On Thu, Jun 6, 2013 at 11:36 AM, William Dunlap  wrote:
> 
>
> Try the following:
>generateABFunction <- function(a, b) {
>   force(a)
>   force(b)
>   function(x) a*x + b
>}
>f12 <- generateABFunction(1, 2)
>f53 <- generateABFunction(5,6)
>f12(10:12) # get 12, 13, 14
>f53(10:12) # get 56, 61, 66
>
> See, e.g., yesterday's discussion under the subject
> "Trying to build up functions with its names by means of lapply"
> on why the force() calls are required.  Read up on R's environments
> to see why f12 and f53 look the same but act differently (hint:
> look at ls.str(environment(f12))).
>
> 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 Liviu Andronic
> > Sent: Thursday, June 06, 2013 8:00 AM
> > To: r-help@r-project.org Help
> > Subject: Re: [R] generate simple function with pre-defined constants
> >
> > On Thu, Jun 6, 2013 at 4:48 PM, Liviu Andronic 
> wrote:
> > > Dear all,
> > > Given:
> > > a <- 2
> > > b <- 3
> > >
> > > I'd like to obtain the following function:
> > > f <- function(x) 2 + 3*x
> > >
> > > but when I do this:
> > > f <- function(x) a + b*x
> > > ##f
> > > ##function(x) a + b*x
> > >
> > > the 'a' and 'b' objects do not get evaluated to their constants. How
> > > could I do that?
> > >
> > I found one solution:
> > a <- 2
> > b <- 3
> > f <- eval(parse(text=paste("function(z)", a, "+ z * ", b)))
> > f
> > ##function(z) 2 + z *  3
> >
> > but I still have nightmares from:
> > > fortune("parse")
> >
> > If the answer is parse() you should usually rethink the question.
> >-- Thomas Lumley
> >   R-help (February 2005)
> >
> > Is there a nicer way to approach this? Thanks,
> > Liviu
> >
> >
> > > Thanks,
> > > Liviu
> > >
> > >
> > > --
> > > Do you know how to read?
> > > http://www.alienetworks.com/srtest.cfm
> > > http://goodies.xfce.org/projects/applications/xfce4-dict#speed-reader
> > > Do you know how to write?
> > > http://garbl.home.comcast.net/~garbl/stylemanual/e.htm#e-mail
> >
> >
> >
> > --
> > Do you know how to read?
> > http://www.alienetworks.com/srtest.cfm
> > http://goodies.xfce.org/projects/applications/xfce4-dict#speed-reader
> > Do you know how to write?
> > http://garbl.home.comcast.net/~garbl/stylemanual/e.htm#e-mail
> >
> > __
> > 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.
>
> ** **
>

[[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] generate simple function with pre-defined constants

2013-06-06 Thread William Dunlap
I said the force was 'required' in the sense that without it
the function will fail to do what you want in some situations.
It doesn't make sense to write a function that you know will
fail sometimes when you know an easy way to make it work
in all situations.

Bill Dunlap
Spotfire, TIBCO Software
wdunlap tibco.com

From: Mark Leeds [mailto:marklee...@gmail.com]
Sent: Thursday, June 06, 2013 8:57 AM
To: William Dunlap
Cc: Liviu Andronic; r-help@r-project.org Help
Subject: Re: [R] generate simple function with pre-defined constants

hi bill: I understand what you're doing but, atleast for this case, I checked 
and you don't need the force this one. it works without it. so, I think the 
force requirement applies only when you're building them up with the lapply. 
but definitely I'm opened to clarification. thanks.



On Thu, Jun 6, 2013 at 11:36 AM, William Dunlap 
mailto:wdun...@tibco.com>> wrote:
Try the following:
   generateABFunction <- function(a, b) {
  force(a)
  force(b)
  function(x) a*x + b
   }
   f12 <- generateABFunction(1, 2)
   f53 <- generateABFunction(5,6)
   f12(10:12) # get 12, 13, 14
   f53(10:12) # get 56, 61, 66

See, e.g., yesterday's discussion under the subject
"Trying to build up functions with its names by means of lapply"
on why the force() calls are required.  Read up on R's environments
to see why f12 and f53 look the same but act differently (hint:
look at ls.str(environment(f12))).

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 Liviu Andronic
> Sent: Thursday, June 06, 2013 8:00 AM
> To: r-help@r-project.org Help
> Subject: Re: [R] generate simple function with pre-defined constants
>
> On Thu, Jun 6, 2013 at 4:48 PM, Liviu Andronic 
> mailto:landronim...@gmail.com>> wrote:
> > Dear all,
> > Given:
> > a <- 2
> > b <- 3
> >
> > I'd like to obtain the following function:
> > f <- function(x) 2 + 3*x
> >
> > but when I do this:
> > f <- function(x) a + b*x
> > ##f
> > ##function(x) a + b*x
> >
> > the 'a' and 'b' objects do not get evaluated to their constants. How
> > could I do that?
> >
> I found one solution:
> a <- 2
> b <- 3
> f <- eval(parse(text=paste("function(z)", a, "+ z * ", b)))
> f
> ##function(z) 2 + z *  3
>
> but I still have nightmares from:
> > fortune("parse")
>
> If the answer is parse() you should usually rethink the question.
>-- Thomas Lumley
>   R-help (February 2005)
>
> Is there a nicer way to approach this? Thanks,
> Liviu
>
>
> > Thanks,
> > Liviu
> >
> >
> > --
> > Do you know how to read?
> > http://www.alienetworks.com/srtest.cfm
> > http://goodies.xfce.org/projects/applications/xfce4-dict#speed-reader
> > Do you know how to write?
> > http://garbl.home.comcast.net/~garbl/stylemanual/e.htm#e-mail
>
>
>
> --
> Do you know how to read?
> http://www.alienetworks.com/srtest.cfm
> http://goodies.xfce.org/projects/applications/xfce4-dict#speed-reader
> Do you know how to write?
> http://garbl.home.comcast.net/~garbl/stylemanual/e.htm#e-mail
>
> __
> 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.


[[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] generate simple function with pre-defined constants

2013-06-06 Thread Mark Leeds
hi bill: I understand what you're doing but, atleast for this case, I
checked and you don't need the force this one. it works without it. so, I
think the force requirement applies only when you're building them up with
the lapply. but definitely I'm opened to clarification. thanks.






On Thu, Jun 6, 2013 at 11:36 AM, William Dunlap  wrote:

> Try the following:
>generateABFunction <- function(a, b) {
>   force(a)
>   force(b)
>   function(x) a*x + b
>}
>f12 <- generateABFunction(1, 2)
>f53 <- generateABFunction(5,6)
>f12(10:12) # get 12, 13, 14
>f53(10:12) # get 56, 61, 66
>
> See, e.g., yesterday's discussion under the subject
> "Trying to build up functions with its names by means of lapply"
> on why the force() calls are required.  Read up on R's environments
> to see why f12 and f53 look the same but act differently (hint:
> look at ls.str(environment(f12))).
>
> 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 Liviu Andronic
> > Sent: Thursday, June 06, 2013 8:00 AM
> > To: r-help@r-project.org Help
> > Subject: Re: [R] generate simple function with pre-defined constants
> >
> > On Thu, Jun 6, 2013 at 4:48 PM, Liviu Andronic 
> wrote:
> > > Dear all,
> > > Given:
> > > a <- 2
> > > b <- 3
> > >
> > > I'd like to obtain the following function:
> > > f <- function(x) 2 + 3*x
> > >
> > > but when I do this:
> > > f <- function(x) a + b*x
> > > ##f
> > > ##function(x) a + b*x
> > >
> > > the 'a' and 'b' objects do not get evaluated to their constants. How
> > > could I do that?
> > >
> > I found one solution:
> > a <- 2
> > b <- 3
> > f <- eval(parse(text=paste("function(z)", a, "+ z * ", b)))
> > f
> > ##function(z) 2 + z *  3
> >
> > but I still have nightmares from:
> > > fortune("parse")
> >
> > If the answer is parse() you should usually rethink the question.
> >-- Thomas Lumley
> >   R-help (February 2005)
> >
> > Is there a nicer way to approach this? Thanks,
> > Liviu
> >
> >
> > > Thanks,
> > > Liviu
> > >
> > >
> > > --
> > > Do you know how to read?
> > > http://www.alienetworks.com/srtest.cfm
> > > http://goodies.xfce.org/projects/applications/xfce4-dict#speed-reader
> > > Do you know how to write?
> > > http://garbl.home.comcast.net/~garbl/stylemanual/e.htm#e-mail
> >
> >
> >
> > --
> > Do you know how to read?
> > http://www.alienetworks.com/srtest.cfm
> > http://goodies.xfce.org/projects/applications/xfce4-dict#speed-reader
> > Do you know how to write?
> > http://garbl.home.comcast.net/~garbl/stylemanual/e.htm#e-mail
> >
> > __
> > 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.
>

[[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] Trying to build up functions with its names by means of lapply

2013-06-06 Thread Bert Gunter
Sergio:

Fair question.

Unfair answer: My personal hangup. To my taste, get() makes R like a
macro language instead of doing functional programming. force() makes
me nervous about how I'm passing arguments. I won't attempt to defend
either of these claims, so feel free to dismiss.

Cheers,
Bert

On Thu, Jun 6, 2013 at 8:28 AM, Julio Sergio  wrote:
> Bert Gunter  gene.com> writes:
>
>> Another equivalent way to do it?
>>
>> f2 <- function(c,nm = "gamma",...)
>> {
>>   probFunc <- paste0(c,nm)
>>   more <- list(...)
>>   function(x)do.call(probFunc,c(x,more))
>> }
>>
>> This avoids the explicit use of get() and force(), I believe, but are
>> there problems here I'm missing?
>>
>
> Thanks Bert. Since I'm relatively new in using this language, a question
> arises in my mind: why should the use of get() and force() be avoided? Is it
> just for syntactic clarity or is there another computational reason?
>
> Best regards,
>
>   -Sergio.
>
> __
> 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.



-- 

Bert Gunter
Genentech Nonclinical Biostatistics

Internal Contact Info:
Phone: 467-7374
Website:
http://pharmadevelopment.roche.com/index/pdb/pdb-functional-groups/pdb-biostatistics/pdb-ncb-home.htm

__
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] Conditional coloring of area between curves

2013-06-06 Thread Roland Pape
Dear list,

I have two time series of temperatures from different sites plotted in the same 
diagram and would like to color the area between the curves differently, 
dependent on whether site 1 is cooler than site 2 (colored let's say in blue) 
or warmer (red). While polygone() works fine to color the enclosed area in 
general, I'm struggling with this conditional coloring. Any help is greatly 
appreciated!

Roland

[[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] generating a bar chart with two axis for co-linear variable

2013-06-06 Thread John Kane
I think we really need to see the code.

John Kane
Kingston ON Canada


> -Original Message-
> From: sudha.krish...@marlabs.com
> Sent: Thu, 6 Jun 2013 06:37:41 +
> To: r-help@r-project.org
> Subject: [R] generating a bar chart with two axis for co-linear variable
> 
> 
> 
> Hello Dimitris,
> 
> 
> 
> I was goggling for some help on Sensitivity vs 1-specificity and saw your
> link.
> 
> 
> 
> I hope you can be of help to me in one of the issue that I am facing in
> generating combo chart(bar chart and plot). I am a novice and have some
> difficulty in getting this logic correct.
> 
> 
> 
> 
> 
> I am give a dataset (I am attaching a sample dataset).
> 
> 
> 
> I am using a barplot() and passing values for percentage frequency and
> the corresponding variables. I am struck here, what my function does is
> only calculate the frequency for the listed variables and not the
> frequency percentage. Is there a method or a script with which I can pass
> the frequency percent and the related values as category columns for x
> axis?
> 
> 
> 
> I will attach the graphs that I have generated so that you can suggest
> the better way.
> 
> 
> 
> Sampledata - Sampledata.txt
> 
> What my function does to calculate the frequency with category names in X
> axis - 1.png
> 
> My requirement is to generate percentage frequency of the variable in y1
> and not the frequency itself. 2.png (where x categories are missing)
> 
> 
> 
> 
> 
> Thanks,
> 
> Sudha Krishnan
> 
> 
> __
> 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.


FREE 3D EARTH SCREENSAVER - Watch the Earth right on your desktop!

__
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] Problem with marginal effects of a multinomial logistic regression

2013-06-06 Thread Leonard Moulin
Hi r-users,

I try to calculate marginal effects of a multinomial logistic regression. To do 
this i use mlogit package and effects() function.

Here is how the procedure works (source : effects() function of mlogit package) 
:
data("Fishing", package = "mlogit")
Fish <- mlogit.data(Fishing, varying = c(2:9), shape = "wide", choice = "mode")
m <- mlogit(mode ~ price | income | catch, data = Fish)
# compute a data.frame containing the mean value of the covariates in
# the sample
z <- with(Fish, data.frame(price = tapply(price, index(m)$alt, mean),
   catch = tapply(catch, index(m)$alt, mean),
   income = mean(income)))
# compute the marginal effects (the second one is an elasticity
effects(m, covariate = "income", data = z)
effects(m, covariate = "price", type = "rr", data = z)
effects(m, covariate = "catch", type = "ar", data = z)
I have no problem with first step (mlogit.data() function). I think my problem 
is on the specification of the multinomial regression.

My regression (for example with three variables) is on the form: Y ~ 0 | X1 + 
X2 + X3. When I try to estimate the marginal effects for model 2 variables no 
problem, however for 3 variables R console returns me the following error: 
"Error in if (rhs% in% c (1, 3)) {: argument is of length zero " (translation 
from error in R console in french).

To understand what is my problem I tried to perform a multinomial regression of 
similar shape on the dataset "Fishing", ie: mode ~ 0 | income + price + catch 
(even if this form has no "economic" sense.) Again the R console returns me the 
same error for 3 variables but manages to estimate these effects for a model 
with two variables ...

This leads me to think that my problem really comes from the specification of 
my multinomial regression ... Do you know how I could find a solution to my 
problem?

Thank you for your help :)

Léonard.
[[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] generate simple function with pre-defined constants

2013-06-06 Thread William Dunlap
Try the following:
   generateABFunction <- function(a, b) {
  force(a)
  force(b)
  function(x) a*x + b
   }
   f12 <- generateABFunction(1, 2)
   f53 <- generateABFunction(5,6)
   f12(10:12) # get 12, 13, 14
   f53(10:12) # get 56, 61, 66

See, e.g., yesterday's discussion under the subject
"Trying to build up functions with its names by means of lapply"
on why the force() calls are required.  Read up on R's environments
to see why f12 and f53 look the same but act differently (hint:
look at ls.str(environment(f12))).

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 Liviu Andronic
> Sent: Thursday, June 06, 2013 8:00 AM
> To: r-help@r-project.org Help
> Subject: Re: [R] generate simple function with pre-defined constants
> 
> On Thu, Jun 6, 2013 at 4:48 PM, Liviu Andronic  wrote:
> > Dear all,
> > Given:
> > a <- 2
> > b <- 3
> >
> > I'd like to obtain the following function:
> > f <- function(x) 2 + 3*x
> >
> > but when I do this:
> > f <- function(x) a + b*x
> > ##f
> > ##function(x) a + b*x
> >
> > the 'a' and 'b' objects do not get evaluated to their constants. How
> > could I do that?
> >
> I found one solution:
> a <- 2
> b <- 3
> f <- eval(parse(text=paste("function(z)", a, "+ z * ", b)))
> f
> ##function(z) 2 + z *  3
> 
> but I still have nightmares from:
> > fortune("parse")
> 
> If the answer is parse() you should usually rethink the question.
>-- Thomas Lumley
>   R-help (February 2005)
> 
> Is there a nicer way to approach this? Thanks,
> Liviu
> 
> 
> > Thanks,
> > Liviu
> >
> >
> > --
> > Do you know how to read?
> > http://www.alienetworks.com/srtest.cfm
> > http://goodies.xfce.org/projects/applications/xfce4-dict#speed-reader
> > Do you know how to write?
> > http://garbl.home.comcast.net/~garbl/stylemanual/e.htm#e-mail
> 
> 
> 
> --
> Do you know how to read?
> http://www.alienetworks.com/srtest.cfm
> http://goodies.xfce.org/projects/applications/xfce4-dict#speed-reader
> Do you know how to write?
> http://garbl.home.comcast.net/~garbl/stylemanual/e.htm#e-mail
> 
> __
> R-help@r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.

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


Re: [R] generate simple function with pre-defined constants

2013-06-06 Thread Pascal Oettli
Hi,

I am not sure what you are looking for. Here are some examples:

foo <- function(a,b,x) a + b*x

> foo
function(a,b,x) a + b*x

a <- 2
b <- 3
x <- 0:10

> foo(a,b,x)
 [1]  2  5  8 11 14 17 20 23 26 29 32

Or

library(polynom)

p1 <- polynomial(c(a,b))
> p1
2 + 3*x

f1 <- as.function(p1)

> f1(x)
 [1]  2  5  8 11 14 17 20 23 26 29 32


Regards,
Pascal



2013/6/6 Liviu Andronic 

> Dear all,
> Given:
> a <- 2
> b <- 3
>
> I'd like to obtain the following function:
> f <- function(x) 2 + 3*x
>
> but when I do this:
> f <- function(x) a + b*x
> ##f
> ##function(x) a + b*x
>
> the 'a' and 'b' objects do not get evaluated to their constants. How
> could I do that?
>
> Thanks,
> Liviu
>
>
> --
> Do you know how to read?
> http://www.alienetworks.com/srtest.cfm
> http://goodies.xfce.org/projects/applications/xfce4-dict#speed-reader
> Do you know how to write?
> http://garbl.home.comcast.net/~garbl/stylemanual/e.htm#e-mail
>
> __
> 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] Trying to build up functions with its names by means of lapply

2013-06-06 Thread Julio Sergio
Bert Gunter  gene.com> writes:

> Another equivalent way to do it?
> 
> f2 <- function(c,nm = "gamma",...)
> {
>   probFunc <- paste0(c,nm)
>   more <- list(...)
>   function(x)do.call(probFunc,c(x,more))
> }
> 
> This avoids the explicit use of get() and force(), I believe, but are
> there problems here I'm missing?
> 

Thanks Bert. Since I'm relatively new in using this language, a question 
arises in my mind: why should the use of get() and force() be avoided? Is it 
just for syntactic clarity or is there another computational reason?

Best regards,

  -Sergio.

__
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] List into table

2013-06-06 Thread David Winsemius

On Jun 5, 2013, at 1:49 AM, Andtrei89 wrote:

> Hi everyone,
> 
> I want to open an XML file in R and it also works.
> But then I get a list of all numbers but it isn't in the type of "numeric".
> 
> I tried this:
> 
> exampleData <-
> "C:/Users/Andreas/Desktop/Thesis/Interzeptionsmodell/Daten/Karlsruhe/Windgeschwindigkeit/1948.xml"
> 
> F <- system.file("exampleData", "1948.xml", package = "XML")
> f <- xmlToDataFrame(exampleData, colClasses = NULL, homogeneous = TRUE,
> collectNames = TRUE,
>stringsAsFactors = default.stringsAsFactors())
> 
> 
> I got this:
> 

You got that when you did … what?

> 1 4.9 5.0 5.3 4.4 2.9 1.9 1.1 1.5 1.6 1.2 0.7 0.8 2.5 2.3 1.6 0.6 1.3 0.8
> 1.3 1.2 3.9 5.0 6.4 6.3 8.4
>   NA  NA  NA  NA  NA  NA  NA  NA  NA  NA  NA  NA  NA  NA  NA  NA  NA  NA 
> NA  NA  NA  NA  NA  NA  NA
> 1 8.4 8.3 8.3 9.1 8.3 7.7 6.8 5.4 4.8 3.6 3.6 3.5 2.6 3.5 3.0 3.1 2.4 3.4
> 3.3 4.8 6.1 5.1 4.2 4.5 4.5
>   NA  NA  NA  NA  NA  NA  NA  NA  NA  NA  NA  NA  NA  NA  NA  NA  NA  NA 
> NA  NA  NA  NA  NA  NA  NA
> 1 5.4 4.7 4.6 3.8 4.0 1.5 0.0 1.1 0.7 0.0 0.0 1.3 1.5 1.6 0.0 0.7 0.6 0.5
> 1.8 1.5 1.6 1.5 1.9 2.5 2.5
> 
> I only want the numbers out of it as numeric ... if I try as.numeric() I got
> a list looking like that:
> 
Probably not a "list". Probably a vector of indexes into a set of factor levels.

> 1 1 1 1 1 1 1 1 1  1  1 1 1 1 1 1 ..
> 
> I would be very thankful for your answer

Hard to know for sure, but try reading the FAQ regarding proper conversion of 
factor varibles to numeric and following its advice.


> --
> View this message in context: 
> http://r.789695.n4.nabble.com/List-into-table-tp4668693.html

Posting from Nabble is likely to cause you problems at some point or another 
due to the tightened spam filters.

> Sent from the R help mailing list archive at Nabble.com.

Nabble is not an archive and it is not the R mailing list.

-- 

David Winsemius
Alameda, CA, USA

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


Re: [R] combining two different matrizes

2013-06-06 Thread arun
Perhaps this also helps:
library(plyr)

 do.call(rbind,alply(aperm(laply(list(A,B),as.matrix),c(1,3,2)),3)) #Using 
Berend's example
   1  2  3  4  5
# [1,] -0.591 -0.934 -0.828  0.012 -0.683
 #[2,]  1.000  6.000 11.000 16.000 21.000
 #[3,]  0.027  1.324 -0.348 -0.223 -0.016
 #[4,]  2.000  7.000 12.000 17.000 22.000
 #[5,] -1.517  0.625 -1.538  0.888 -0.443
 #[6,]  3.000  8.000 13.000 18.000 23.000
 #[7,] -1.363 -0.046 -0.256 -0.592  0.353
 #[8,]  4.000  9.000 14.000 19.000 24.000
 #[9,]  1.178 -1.004 -1.150 -0.656  0.073
#[10,]  5.000 10.000 15.000 20.000 25.000


A.K.

- Original Message -
From: Berend Hasselman 
To: ThomasH 
Cc: r-help@r-project.org
Sent: Thursday, June 6, 2013 10:04 AM
Subject: Re: [R] combining two different matrizes


On 05-06-2013, at 23:56, ThomasH  wrote:

> 
> Hello together,
> 
> this is ma first post, so please aplogize me if post this in the wrong
> section.
> 
> I have problem concerning ma two matrizes.
> 
> After a regressione and so on, I got two matrizes
> 
> Matrixres contains the results of ma calculation.
> 
> Matrixr contains my detiene, which where Aldo used for the regression.
> 
> Please ser the following code:
> 
> #Datei einlesen
> residual = read.csv2("E:***Input-R_Renditen.csv",header=TRUE, sep=";")
> 
> 
> #Aktientitel
> alist <- list()
> for (a in 2:11){
> 
> 
> #Regression
>    #Länge Gesamtzeit
>    t <- 243
>    tx <- t-59
> 
>    #Länge Regression
>    reglist <- list()
>    for (i in 1:tx){
>    j <- i+59
> 
>    #RegressionsVariable
>    x = residual[i:j,a]
>    rm = residual[i:j,12]
>    smb = residual[i:j,13]
>    hml = residual[i:j,14]
>    rf = residual[i:j,15]
> 
>    #Überschussrenditen
>    ex=x-rf
>    erm=rm-rf
> 
>    #Regression
>    reg <- lm(ex~erm+smb+hml)
>    reglist[[i]] <- coef(reg)
>    
> 
> #Berechnung Residuum
>       #Residual Berechnung
>       rx = residual[(j-5):j,a]
>       rrm = residual[(j-5):j,12]
>       rsmb = residual[(j-5):j,13]
>       rhml = residual[(j-5):j,14]
>       rrf = residual[(j-5):j,15]
> 
>       rex = rx-rrf
>       rerm = rrm-rrf
> 
>       #Berechnung
>       res <-
> sum(rex-(reglist[[i]][2]*rerm+reglist[[i]][3]*rsmb+reglist[[i]][4]*rhml))/sd(rex-(reglist[[i]][2]*rerm+reglist[[i]][3]*rsmb+reglist[[i]][4]*rhml))
>       reglist[[i]] <- res
> }  
>    
> 
> #Residuen auf alle Aktien
>    alist[[a]] <- reglist
> }
>    
> #Matrix mit Residuen
>    matrixres <- do.call(cbind,alist)
> 
> #Spaltennamen/Zeilennamen
> s<- names(residual)[2:11]
> colnames(matrixres)<-s
> 
> #RenditeMatrix
> matrixr <- do.call(cbind,residual[60:243,2:11])
> 
> 
> Now I want to combines  the two matrizes in the following way:
> 
> Under every row of matrixres should stand the row of matrixr for excample:
> 
> Matrixres row1
> Matrixr row1
> Matrixres row2
> Matrixr row 2
> 
> Can anybody help me? I was working on this problem the whole day, but have
> no idea.

Something like this (assuming A and B have the same dimensions)

set.seed(11)
A <- matrix(round(rnorm(25),3),nrow=5)
B <- matrix(1:25,nrow=5)

C <- matrix(0,nrow=2*nrow(A),ncol=ncol(A))

crows <- seq.int(from=1,to=2*nrow(A),by=2)
C[crows,] <- A
C[crows+1,] <- B

Some friendly advice: get someone to check your English before sending a mail 
to the list.

Berend

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

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


Re: [R] SPlus script

2013-06-06 Thread Scott Raynaud
I actually had tried placing arguments in the call but it didn't work.   
However, I did 
not think about writing it to a variable and printing.  That seems to have done 
the 
trick.  Funny, I don't remember having to do that before, but that's not 
surprising.

Anyway, thanks for helping to diagnose and fix the problem.

- Original Message -
From: Ista Zahn 
To: Scott Raynaud 
Cc: William Dunlap ; "r-help@r-project.org" 

Sent: Thursday, June 6, 2013 9:15 AM
Subject: Re: [R] SPlus script

Presumably something like

r <- sshc(50)
print(r)

But if you were getting output before than you already have a script
that does something like this. It would be better to find it...

Best,
Ista

On Thu, Jun 6, 2013 at 9:02 AM, Scott Raynaud  wrote:
> Ok.  Now I see that the sshc function is not being called.  Thanks for 
> pointing that out.
> I'm not certain about the solution, however.  I tried putting call("sshc") at 
> the end of the
> program, but nothing happened.  My memory about all of this is fuzzy.  
> Suggestions
> on how to call the function appreciated.
>
> - Original Message -
> From: William Dunlap 
> To: Scott Raynaud ; "r-help@r-project.org" 
> 
> Cc:
> Sent: Wednesday, June 5, 2013 2:17 PM
> Subject: RE: [R] SPlus script
>
> Both the R and S+ versions (which seem to differ only in the use of _ for 
> assignment
> in the S+ version) do nothing but define some functions.  You would not 
> expect any
> printed output unless you used those functions on some data.  Is there 
> another script
> that does that?
>
> 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 Scott Raynaud
>> Sent: Wednesday, June 05, 2013 6:21 AM
>> To: r-help@r-project.org
>> Subject: [R] SPlus script
>>
>> This originally was an SPlus script that I modifeid about a year-and-a-half 
>> ago.  It worked
>> perfectly then.  Now I can't get any output despite not receiving an error 
>> message.  I'm
>> providing the SPLUS script as a reference.  I'm running R15.2.2.  Any help 
>> appreciated.
>>
>> MY
>> MODIFICATION***
>> **
>> ## sshc.ssc: sample size calculation for historical control studies
>> ## J. Jack Lee (jj...@mdanderson.org) and Chi-hong Tseng
>> ## Department of Biostatistics, Univ. of Texas M.D. Anderson Cancer Center
>> ##
>> ## 3/1/99
>> ## updated 6/7/00: add loess
>> ##--
>>  Required Input:
>> #
>> # rc    number of response in historical control group
>> # nc    sample size in historical control
>> # d      target improvement = Pe - Pc
>> # method 1=method based on the randomized design
>> #        2=Makuch & Simon method (Makuch RW, Simon RM. Sample size 
>> considerations
>> #          for non-randomized comparative studies. J of Chron Dis 1980; 
>> 3:175-181.
>> #        3=uniform power method
>>  optional Input:
>> #
>> # alpha  size of the test
>> # power  desired power of the test
>> # tol    convergence criterion for methods 1 & 2 in terms of sample size
>> # tol1  convergence criterion for method 3 at any given obs Rc in terms of 
>> difference
>> #          of expected power from target
>> # tol2  overall convergence criterion for method 3 as the max absolute 
>> deviation
>> #          of expected power from target for all Rc
>> # cc    range of multiplicative constant applied to the initial values ne
>> # l.span smoothing constant for loess
>> #
>> # Note:  rc is required for methods 1 and 2 but not 3
>> #        method 3 return the sample size need for rc=0 to (1-d)*nc
>> #
>>  Output
>> # for methdos 1 & 2: return the sample size needed for the experimental 
>> group (1
>> number)
>> #                    for given rc, nc, d, alpha, and power
>> # for method 3:      return the profile of sample size needed for given nc, 
>> d, alpha, and
>> power
>> #                    vector $ne contains the sample size corresponding to 
>> rc=0, 1, 2, ... nc*(1-d)
>> #                    vector $Ep contains the expected power corresponding to
>> #                      the true pc = (0, 1, 2, ..., nc*(1-d)) / nc
>> #
>> #--
>> sshc<-function(rc, nc=1092, d=.085779816, method=3, alpha=0.05, power=0.8,
>>              tol=0.01, tol1=.0001, tol2=.005, cc=c(.1,2), l.span=.5)
>> {
>> ### for method 1
>> if (method==1) {
>>  ne1<-ss.rand(rc,nc,d,alpha=.05,power=.8,tol=.01)
>>  return(ne=ne1)
>>                }
>> ### for method 2
>> if (method==2) {
>> ne<-nc
>> ne1<-nc+50
>> while(abs(ne-ne1)>tol & ne1<10){
>> ne<-ne1
>> pe<-d+rc/nc
>> ne1<-nef(rc,nc,pe*ne,ne,alpha,power)
>> ## if(is.na(ne1)) print(paste('rc=',rc,',nc=',nc,',pe=',pe,',ne=',ne))
>> }
>> if (ne1>10) return(NA)
>> else return(ne=ne1)
>> }
>> ### for method

Re: [R] combining two different matrizes

2013-06-06 Thread Berend Hasselman

On 06-06-2013, at 17:04, David Carlson  wrote:

> You didn't give us data, but this may give you enough to solve your problem:
> 
>> set.seed(42)
>> nrows <- 6
>> ncols <- 5
>> mat1 <- matrix(sample.int(100, 30), nrows, ncols)
>> mat1
> …..
>> newmat <- matrix(rbind(as.vector(mat1), as.vector(mat2)), nrows*2, ncols)

Nice.

Berend

__
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] generate simple function with pre-defined constants

2013-06-06 Thread Liviu Andronic
On Thu, Jun 6, 2013 at 5:03 PM, arun  wrote:
> HI,
> Not sure I understand your question:
>  a <- 2
>  b <- 3
>  f1<- function(x) a+b*x
>
I don't want the function to depend on the objects a and b, but
instead use the values of those objects (I do this within a function).

Liviu


>  f1(2)
> #[1] 8
>  f1(3)
> #[1] 11
>  f<- function(x) 2+3*x
>  f(2)
> #[1] 8
>  f(3)
> #[1] 11
>
>
> A.K.
>
>   sessionInfo()
> R version 3.0.0 (2013-04-03)
> Platform: x86_64-unknown-linux-gnu (64-bit)
>
> locale:
>  [1] LC_CTYPE=en_CA.UTF-8   LC_NUMERIC=C
>  [3] LC_TIME=en_CA.UTF-8LC_COLLATE=en_CA.UTF-8
>  [5] LC_MONETARY=en_CA.UTF-8LC_MESSAGES=en_CA.UTF-8
>  [7] LC_PAPER=C LC_NAME=C
>  [9] LC_ADDRESS=C   LC_TELEPHONE=C
> [11] LC_MEASUREMENT=en_CA.UTF-8 LC_IDENTIFICATION=C
>
> attached base packages:
> [1] stats graphics  grDevices utils datasets  methods   base
>
> other attached packages:
> [1] arrayhelpers_0.76-20120816 abind_1.4-0
> [3] plyr_1.8   stringr_0.6.2
> [5] reshape2_1.2.2
>
> loaded via a namespace (and not attached):
> [1] tools_3.0.0
>
>
> - Original Message -
> From: Liviu Andronic 
> To: "r-help@r-project.org Help" 
> Cc:
> Sent: Thursday, June 6, 2013 10:48 AM
> Subject: [R] generate simple function with pre-defined constants
>
> Dear all,
> Given:
> a <- 2
> b <- 3
>
> I'd like to obtain the following function:
> f <- function(x) 2 + 3*x
>
> but when I do this:
> f <- function(x) a + b*x
> ##f
> ##function(x) a + b*x
>
> the 'a' and 'b' objects do not get evaluated to their constants. How
> could I do that?
>
> Thanks,
> Liviu
>
>
> --
> Do you know how to read?
> http://www.alienetworks.com/srtest.cfm
> http://goodies.xfce.org/projects/applications/xfce4-dict#speed-reader
> Do you know how to write?
> http://garbl.home.comcast.net/~garbl/stylemanual/e.htm#e-mail
>
> __
> 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.
>



-- 
Do you know how to read?
http://www.alienetworks.com/srtest.cfm
http://goodies.xfce.org/projects/applications/xfce4-dict#speed-reader
Do you know how to write?
http://garbl.home.comcast.net/~garbl/stylemanual/e.htm#e-mail

__
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] Not sure this is something R could do but it feels like it should be.

2013-06-06 Thread Marc Schwartz

On Jun 6, 2013, at 10:03 AM, Polwart Calum (COUNTY DURHAM AND DARLINGTON NHS 
FOUNDATION TRUST)  wrote:

> Some colleagues nationally have developed a system which means they can pick 
> the optimal sets of doses for a drug.  The system could apply to a number of 
> drugs.  But the actual doses might vary.  To try and explain this in terms 
> that the average Joe on the street might understand if you have some 
> amoxicillin antibiotic for a chest infection the normal dose for an adult is 
> 250 to 500mg increased to maybe 1000mg in severe cases.
> 
> For a child it is dosed from a liquid and people usually go from 62.5mg, 
> 125mg to 250mg although you could measure any volume you wanted.
> 
> What this new method has developed is a means to pick the "right" standard 
> doses so what above is 62.5, 125, 250, 500, 1000.  However the method they've 
> used is really engineered about ensure the jump between doses is correct - 
> you'll notice that the list above is a doubling up method.
> 
> But you can also have a doubling up method that went 50, 100, 200, 400, 800, 
> 1600  and pretty much as many as you can think of depending on your starting 
> point and there is no scientific means to pick that starting point.  So 
> colleagues have developed their rather more complex equivalent of the 
> doubling method to determine the doses they need but they need to know if 
> they should start at 40, 50, 62.5 or some other number.
> 
> Once they have the starting number they can calculate all the other doses.  I 
> realise R can do that, and I realise using a loop of possible starting 
> numbers it can build all those options.
> 
> Each patient then has a theoretical dose they should get lets say that's 
> 10mg/kg and you might treat patients from 5 to 120kg.  They are then looking 
> to calculate the variance for each dose range so if we take the 50, 100, 200, 
> 400 model and said you'd give 50mg to anyone needing 0?? to 75mg 100mg to 
> anyone needing 76 - 150mg etc... from there they are taking that range and 
> saying that's a 31% overdose to a 33% underdose.  Then they want to find if 
> there is a starting number which minimises the extent of under and 
> overdosing...
> 
> Anyone know of an existing stats function in R that can easily do that and 
> almost then report from some inputs a single number that is the "best fit"?
> 
> Calum


The first place I would start is with the two relevant CRAN Task Views:

  http://cran.r-project.org/web/views/ClinicalTrials.html

and

  http://cran.r-project.org/web/views/Pharmacokinetics.html


There is also another package not listed above that might be relevant:

  http://cran.r-project.org/web/packages/scaRabee/


Regards,

Marc Schwartz

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


Re: [R] generate simple function with pre-defined constants

2013-06-06 Thread arun
HI,
Not sure I understand your question:
 a <- 2
 b <- 3
 f1<- function(x) a+b*x
 f1(2)
#[1] 8
 f1(3)
#[1] 11
 f<- function(x) 2+3*x
 f(2)
#[1] 8
 f(3)
#[1] 11


A.K.

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

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

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

other attached packages:
[1] arrayhelpers_0.76-20120816 abind_1.4-0   
[3] plyr_1.8   stringr_0.6.2 
[5] reshape2_1.2.2    

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


- Original Message -
From: Liviu Andronic 
To: "r-help@r-project.org Help" 
Cc: 
Sent: Thursday, June 6, 2013 10:48 AM
Subject: [R] generate simple function with pre-defined constants

Dear all,
Given:
a <- 2
b <- 3

I'd like to obtain the following function:
f <- function(x) 2 + 3*x

but when I do this:
f <- function(x) a + b*x
##f
##function(x) a + b*x

the 'a' and 'b' objects do not get evaluated to their constants. How
could I do that?

Thanks,
Liviu


-- 
Do you know how to read?
http://www.alienetworks.com/srtest.cfm
http://goodies.xfce.org/projects/applications/xfce4-dict#speed-reader
Do you know how to write?
http://garbl.home.comcast.net/~garbl/stylemanual/e.htm#e-mail

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


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


Re: [R] generate simple function with pre-defined constants

2013-06-06 Thread Liviu Andronic
On Thu, Jun 6, 2013 at 4:48 PM, Liviu Andronic  wrote:
> Dear all,
> Given:
> a <- 2
> b <- 3
>
> I'd like to obtain the following function:
> f <- function(x) 2 + 3*x
>
> but when I do this:
> f <- function(x) a + b*x
> ##f
> ##function(x) a + b*x
>
> the 'a' and 'b' objects do not get evaluated to their constants. How
> could I do that?
>
I found one solution:
a <- 2
b <- 3
f <- eval(parse(text=paste("function(z)", a, "+ z * ", b)))
f
##function(z) 2 + z *  3

but I still have nightmares from:
> fortune("parse")

If the answer is parse() you should usually rethink the question.
   -- Thomas Lumley
  R-help (February 2005)

Is there a nicer way to approach this? Thanks,
Liviu


> Thanks,
> Liviu
>
>
> --
> Do you know how to read?
> http://www.alienetworks.com/srtest.cfm
> http://goodies.xfce.org/projects/applications/xfce4-dict#speed-reader
> Do you know how to write?
> http://garbl.home.comcast.net/~garbl/stylemanual/e.htm#e-mail



-- 
Do you know how to read?
http://www.alienetworks.com/srtest.cfm
http://goodies.xfce.org/projects/applications/xfce4-dict#speed-reader
Do you know how to write?
http://garbl.home.comcast.net/~garbl/stylemanual/e.htm#e-mail

__
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] combining two different matrizes

2013-06-06 Thread David Carlson
You didn't give us data, but this may give you enough to solve your problem:

> set.seed(42)
> nrows <- 6
> ncols <- 5
> mat1 <- matrix(sample.int(100, 30), nrows, ncols)
> mat1
 [,1] [,2] [,3] [,4] [,5]
[1,]   92   70   83   397
[2,]   93   13   23   46   82
[3,]   29   61   40   73   98
[4,]   81   65   80   11   67
[5,]   62   42   88   78   33
[6,]   50   91   10   85   60
> mat2 <- round(prop.table(mat1)*100, 2)
> mat2
 [,1] [,2] [,3] [,4] [,5]
[1,] 5.25 4.00 4.74 2.23 0.40
[2,] 5.31 0.74 1.31 2.63 4.68
[3,] 1.66 3.48 2.28 4.17 5.59
[4,] 4.62 3.71 4.57 0.63 3.82
[5,] 3.54 2.40 5.02 4.45 1.88
[6,] 2.85 5.19 0.57 4.85 3.42
> newmat <- matrix(rbind(as.vector(mat1), as.vector(mat2)), nrows*2, ncols)
> newmat
   [,1]  [,2]  [,3]  [,4]  [,5]
 [1,] 92.00 70.00 83.00 39.00  7.00
 [2,]  5.25  4.00  4.74  2.23  0.40
 [3,] 93.00 13.00 23.00 46.00 82.00
 [4,]  5.31  0.74  1.31  2.63  4.68
 [5,] 29.00 61.00 40.00 73.00 98.00
 [6,]  1.66  3.48  2.28  4.17  5.59
 [7,] 81.00 65.00 80.00 11.00 67.00
 [8,]  4.62  3.71  4.57  0.63  3.82
 [9,] 62.00 42.00 88.00 78.00 33.00
[10,]  3.54  2.40  5.02  4.45  1.88
[11,] 50.00 91.00 10.00 85.00 60.00
[12,]  2.85  5.19  0.57  4.85  3.42

-
David L Carlson
Associate Professor of Anthropology
Texas A&M University
College Station, TX 77840-4352

-Original Message-
From: r-help-boun...@r-project.org [mailto:r-help-boun...@r-project.org] On 
Behalf Of ThomasH
Sent: Wednesday, June 5, 2013 4:56 PM
To: r-help@r-project.org
Subject: [R] combining two different matrizes


Hello together,

this is ma first post, so please aplogize me if post this in the wrong
section.

I have problem concerning ma two matrizes.

After a regressione and so on, I got two matrizes

Matrixres contains the results of ma calculation.

Matrixr contains my detiene, which where Aldo used for the regression.

Please ser the following code:

#Datei einlesen
residual = read.csv2("E:***Input-R_Renditen.csv",header=TRUE, sep=";")


#Aktientitel
alist <- list()
for (a in 2:11){


#Regression
   #Länge Gesamtzeit
   t <- 243
   tx <- t-59

   #Länge Regression
   reglist <- list()
   for (i in 1:tx){
   j <- i+59

   #RegressionsVariable
   x = residual[i:j,a]
   rm = residual[i:j,12]
   smb = residual[i:j,13]
   hml = residual[i:j,14]
   rf = residual[i:j,15]

   #Überschussrenditen
   ex=x-rf
   erm=rm-rf

   #Regression
   reg <- lm(ex~erm+smb+hml)
   reglist[[i]] <- coef(reg)
   

#Berechnung Residuum
  #Residual Berechnung
  rx = residual[(j-5):j,a]
  rrm = residual[(j-5):j,12]
  rsmb = residual[(j-5):j,13]
  rhml = residual[(j-5):j,14]
  rrf = residual[(j-5):j,15]

  rex = rx-rrf
  rerm = rrm-rrf

  #Berechnung
  res <-
sum(rex-(reglist[[i]][2]*rerm+reglist[[i]][3]*rsmb+reglist[[i]][4]*rhml))/sd(rex-(reglist[[i]][2]*rerm+reglist[[i]][3]*rsmb+reglist[[i]][4]*rhml))
  reglist[[i]] <- res
}   
   

#Residuen auf alle Aktien
   alist[[a]] <- reglist
}
   
#Matrix mit Residuen
   matrixres <- do.call(cbind,alist)

#Spaltennamen/Zeilennamen
s<- names(residual)[2:11]
colnames(matrixres)<-s

#RenditeMatrix
matrixr <- do.call(cbind,residual[60:243,2:11])


Now I want to combines  the two matrizes in the following way:

Under every row of matrixres should stand the row of matrixr for excample:

Matrixres row1
Matrixr row1
Matrixres row2
Matrixr row 2

Can anybody help me? I was working on this problem the whole day, but have
no idea.

Thanks alot
Thomash




--
View this message in context: 
http://r.789695.n4.nabble.com/combining-two-different-matrizes-tp4668766.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.


[R] Not sure this is something R could do but it feels like it should be.

2013-06-06 Thread Polwart Calum (COUNTY DURHAM AND DARLINGTON NHS FOUNDATION TRUST)
Some colleagues nationally have developed a system which means they can pick 
the optimal sets of doses for a drug.  The system could apply to a number of 
drugs.  But the actual doses might vary.  To try and explain this in terms that 
the average Joe on the street might understand if you have some amoxicillin 
antibiotic for a chest infection the normal dose for an adult is 250 to 500mg 
increased to maybe 1000mg in severe cases.

For a child it is dosed from a liquid and people usually go from 62.5mg, 125mg 
to 250mg although you could measure any volume you wanted.

What this new method has developed is a means to pick the "right" standard 
doses so what above is 62.5, 125, 250, 500, 1000.  However the method they've 
used is really engineered about ensure the jump between doses is correct - 
you'll notice that the list above is a doubling up method.

But you can also have a doubling up method that went 50, 100, 200, 400, 800, 
1600  and pretty much as many as you can think of depending on your starting 
point and there is no scientific means to pick that starting point.  So 
colleagues have developed their rather more complex equivalent of the doubling 
method to determine the doses they need but they need to know if they should 
start at 40, 50, 62.5 or some other number.

Once they have the starting number they can calculate all the other doses.  I 
realise R can do that, and I realise using a loop of possible starting numbers 
it can build all those options.

Each patient then has a theoretical dose they should get lets say that's 
10mg/kg and you might treat patients from 5 to 120kg.  They are then looking to 
calculate the variance for each dose range so if we take the 50, 100, 200, 400 
model and said you'd give 50mg to anyone needing 0?? to 75mg 100mg to anyone 
needing 76 - 150mg etc... from there they are taking that range and saying 
that's a 31% overdose to a 33% underdose.  Then they want to find if there is a 
starting number which minimises the extent of under and overdosing...

Anyone know of an existing stats function in R that can easily do that and 
almost then report from some inputs a single number that is the "best fit"?

Calum



This message may contain confidential information. If yo...{{dropped:22}}

__
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] rJava is not loading

2013-06-06 Thread Prof Brian Ripley

On 06/06/2013 15:52, Dimitri Liakhovitski wrote:

Thank you very much, Brian.
It's clearly christal clear - but one needs a christal ball to realize
that! :-)
So I learned: architecture = bitness


On current Windows, yes.
But not on OS X nor Linux nor Solaris nor FreeBSD 

It indicates the type of CPU, or more precisely the instruction set the 
CPU is running.  x86_64 CPUs (as used by most current Windows boxes) can 
emulate i386 CPUs.  So the architectures for R on Windows are


i386 (32-bit)
x64 (64-bit, the Windows name for x86_64).


Dimitri


On Thu, Jun 6, 2013 at 2:19 AM, Prof Brian Ripley mailto:rip...@stats.ox.ac.uk>> wrote:

On 06/06/2013 00:38, Dimitri Liakhovitski wrote:

Hello!
I installed rJava and am trying to load it.

library(rJava)

Error : .onLoad failed in loadNamespace() for 'rJava', details:
call: fun(libname, pkgname)
error: No CurrentVersion entry in Software/JavaSoft
registry! Try
re-installing Java and make sure R and Java have matching
architectures.
Error: package or namespace load failed for ‘rJava’

Any idea why?
Background info:


But the posting guide asked for the output from sessionInfo(): we
want to know what R knows it is running as, not why you think.


Windows 7, 64-bit
R version 3.0.1 (for 64-bit)
I just installed the lastest Java: Java 7 Update 21


For the right architecture?

The message is crystal clear: you do not have a Java installed
matching your R.  Most likely your Java is 32-bit and your R 64-bit
or v.v.



--
Brian D. Ripley, rip...@stats.ox.ac.uk 
Professor of Applied Statistics,
http://www.stats.ox.ac.uk/~__ripley/

University of Oxford, Tel: +44 1865 272861
 (self)
1 South Parks Road, +44 1865 272866  (PA)
Oxford OX1 3TG, UKFax: +44 1865 272595





--
Dimitri Liakhovitski



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

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


Re: [R] rJava is not loading

2013-06-06 Thread Dimitri Liakhovitski
Thank you very much, Brian.

It's clearly christal clear - but one needs a christal ball to realize
that! :-)

So I learned: architecture = bitness

Dimitri



On Thu, Jun 6, 2013 at 2:19 AM, Prof Brian Ripley wrote:

> On 06/06/2013 00:38, Dimitri Liakhovitski wrote:
>
>> Hello!
>> I installed rJava and am trying to load it.
>>
>> library(rJava)
>>
>> Error : .onLoad failed in loadNamespace() for 'rJava', details:
>>call: fun(libname, pkgname)
>>error: No CurrentVersion entry in Software/JavaSoft registry! Try
>> re-installing Java and make sure R and Java have matching architectures.
>> Error: package or namespace load failed for ‘rJava’
>>
>> Any idea why?
>> Background info:
>>
>
> But the posting guide asked for the output from sessionInfo(): we want to
> know what R knows it is running as, not why you think.
>
>
>  Windows 7, 64-bit
>> R version 3.0.1 (for 64-bit)
>> I just installed the lastest Java: Java 7 Update 21
>>
>
> For the right architecture?
>
> The message is crystal clear: you do not have a Java installed matching
> your R.  Most likely your Java is 32-bit and your R 64-bit or v.v.
>
>
>
> --
> Brian D. Ripley,  rip...@stats.ox.ac.uk
> Professor of Applied Statistics,  
> http://www.stats.ox.ac.uk/~**ripley/
> University of Oxford, Tel:  +44 1865 272861 (self)
> 1 South Parks Road, +44 1865 272866 (PA)
> Oxford OX1 3TG, UKFax:  +44 1865 272595
>



-- 
Dimitri Liakhovitski

[[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] generate simple function with pre-defined constants

2013-06-06 Thread Liviu Andronic
Dear all,
Given:
a <- 2
b <- 3

I'd like to obtain the following function:
f <- function(x) 2 + 3*x

but when I do this:
f <- function(x) a + b*x
##f
##function(x) a + b*x

the 'a' and 'b' objects do not get evaluated to their constants. How
could I do that?

Thanks,
Liviu


-- 
Do you know how to read?
http://www.alienetworks.com/srtest.cfm
http://goodies.xfce.org/projects/applications/xfce4-dict#speed-reader
Do you know how to write?
http://garbl.home.comcast.net/~garbl/stylemanual/e.htm#e-mail

__
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] SPlus script

2013-06-06 Thread Ista Zahn
Presumably something like

r <- sshc(50)
print(r)

But if you were getting output before than you already have a script
that does something like this. It would be better to find it...

Best,
Ista

On Thu, Jun 6, 2013 at 9:02 AM, Scott Raynaud  wrote:
> Ok.  Now I see that the sshc function is not being called.  Thanks for 
> pointing that out.
> I'm not certain about the solution, however.  I tried putting call("sshc") at 
> the end of the
> program, but nothing happened.  My memory about all of this is fuzzy.  
> Suggestions
> on how to call the function appreciated.
>
> - Original Message -
> From: William Dunlap 
> To: Scott Raynaud ; "r-help@r-project.org" 
> 
> Cc:
> Sent: Wednesday, June 5, 2013 2:17 PM
> Subject: RE: [R] SPlus script
>
> Both the R and S+ versions (which seem to differ only in the use of _ for 
> assignment
> in the S+ version) do nothing but define some functions.  You would not 
> expect any
> printed output unless you used those functions on some data.  Is there 
> another script
> that does that?
>
> 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 Scott Raynaud
>> Sent: Wednesday, June 05, 2013 6:21 AM
>> To: r-help@r-project.org
>> Subject: [R] SPlus script
>>
>> This originally was an SPlus script that I modifeid about a year-and-a-half 
>> ago.  It worked
>> perfectly then.  Now I can't get any output despite not receiving an error 
>> message.  I'm
>> providing the SPLUS script as a reference.  I'm running R15.2.2.  Any help 
>> appreciated.
>>
>> MY
>> MODIFICATION***
>> **
>> ## sshc.ssc: sample size calculation for historical control studies
>> ## J. Jack Lee (jj...@mdanderson.org) and Chi-hong Tseng
>> ## Department of Biostatistics, Univ. of Texas M.D. Anderson Cancer Center
>> ##
>> ## 3/1/99
>> ## updated 6/7/00: add loess
>> ##--
>>  Required Input:
>> #
>> # rc number of response in historical control group
>> # nc sample size in historical control
>> # d  target improvement = Pe - Pc
>> # method 1=method based on the randomized design
>> #2=Makuch & Simon method (Makuch RW, Simon RM. Sample size 
>> considerations
>> #  for non-randomized comparative studies. J of Chron Dis 1980; 
>> 3:175-181.
>> #3=uniform power method
>>  optional Input:
>> #
>> # alpha  size of the test
>> # power  desired power of the test
>> # tolconvergence criterion for methods 1 & 2 in terms of sample size
>> # tol1   convergence criterion for method 3 at any given obs Rc in terms of 
>> difference
>> #  of expected power from target
>> # tol2   overall convergence criterion for method 3 as the max absolute 
>> deviation
>> #  of expected power from target for all Rc
>> # cc range of multiplicative constant applied to the initial values ne
>> # l.span smoothing constant for loess
>> #
>> # Note:  rc is required for methods 1 and 2 but not 3
>> #method 3 return the sample size need for rc=0 to (1-d)*nc
>> #
>>  Output
>> # for methdos 1 & 2: return the sample size needed for the experimental 
>> group (1
>> number)
>> #for given rc, nc, d, alpha, and power
>> # for method 3:  return the profile of sample size needed for given nc, 
>> d, alpha, and
>> power
>> #vector $ne contains the sample size corresponding to 
>> rc=0, 1, 2, ... nc*(1-d)
>> #vector $Ep contains the expected power corresponding to
>> #  the true pc = (0, 1, 2, ..., nc*(1-d)) / nc
>> #
>> #--
>> sshc<-function(rc, nc=1092, d=.085779816, method=3, alpha=0.05, power=0.8,
>>   tol=0.01, tol1=.0001, tol2=.005, cc=c(.1,2), l.span=.5)
>> {
>> ### for method 1
>> if (method==1) {
>>  ne1<-ss.rand(rc,nc,d,alpha=.05,power=.8,tol=.01)
>>  return(ne=ne1)
>>}
>> ### for method 2
>> if (method==2) {
>> ne<-nc
>> ne1<-nc+50
>> while(abs(ne-ne1)>tol & ne1<10){
>> ne<-ne1
>> pe<-d+rc/nc
>> ne1<-nef(rc,nc,pe*ne,ne,alpha,power)
>> ## if(is.na(ne1)) print(paste('rc=',rc,',nc=',nc,',pe=',pe,',ne=',ne))
>> }
>> if (ne1>10) return(NA)
>> else return(ne=ne1)
>> }
>> ### for method 3
>> if (method==3) {
>> if (tol1 > tol2/10) tol1<-tol2/10
>> ncstar<-(1-d)*nc
>> pc<-(0:ncstar)/nc
>> ne<-rep(NA,ncstar + 1)
>> for (i in (0:ncstar))
>> { ne[i+1]<-ss.rand(i,nc,d,alpha=.05,power=.8,tol=.01)
>> }
>> plot(pc,ne,type='l',ylim=c(0,max(ne)*1.5))
>> ans<-c.searchd(nc, d, ne, alpha, power, cc, tol1)
>> ### check overall absolute deviance
>> old.abs.dev<-sum(abs(ans$Ep-power))
>> ##bad<-0
>> print(round(ans$Ep,4))
>> print(round(ans$ne,2))
>> lines(pc,ans$ne,lty=1,col=8)
>> old.ne<

[R] Autocorrelation and normal distribution of gaps for ping requests in an unstable network

2013-06-06 Thread Ramon Hofer
Hi all

I have a powerline network connection which I'm investigating.
The test network contains some nodes to which I ping from one host.
The source host is always the same and I split the data to get files
for each connection.
A lot of ping requests get lost and I'm trying to plot an
autocorrelation of the data.

Here's an example log:
http://people.ee.ethz.ch/~hoferr/download/data-20130603-192.168.72.33.csv

I tried to plot the autocorrelation graph:
 acf(A$pingRTT.ms.)
which didn't work because of missing ping values. I found a post at
stackoverflow [1] where they suggest to use coredata which didn't work
for me. They also suggest to use "na.action = na.omit" or "na.action =
na.pass". The second option works for me.

With these two commands I can draw an autocorrelation graph.
 A <- read.csv('data-20130603-192.168.72.33.csv')
 acf(A$pingRTT.ms., na.action = na.pass)

But they also warn that:
"acf works on regularly spaced data so acf first expands the time
series to a regularly spaced one inserting NAs as needed to make it
regularly spaced."
This seems to me as if it introduces new periods of time where there's
no ping value and thus no connection which means the autocorrelation
graph I get is nonsense.
Is my fear for no reason or is there a way to get a meaningful plot?


I'd also like to plot a histogram with normal curve like the example
from statmethods [2].
In their example they have the data directly available.
In my case I need to prepare my data to get a list of gaps. E.g.

 TimestampStart,GapLength
 2013-06-03_15:20:25.374096766,16.2s
 2013-06-03_15:22:13.944293504,37.5s
 ...

My plan is to program a loop like

 A$Timestamp <- strptime(as.character(A$Timestamp), "%Y-%m-%d_%H:%M:%S")
 B <- matrix(nrow = 0, ncol = 2)
 colnames(B) <- c("TimestampStart","GapLength[s]")
 j <- 1
 gap.start <- A$Timestamp[0]
 for(i in 2:length(A$Timestamp)) 
 { #For all rows
  if(is.na(A$pingRTT.ms.[i]))
  { #Currently no connection
   if(!is.na(A$pingRTT.ms.[i-1]))
   { #Connection lost now
gap.start <- i
   } 
   else if(!is.na(A$pingRTT.ms.[i+1])) 
   { # Connection restores next time
gap.end <- i+1
B <-
 rbind( B,
  c(
   A$Timestamp[gap.start],
   A$Timestamp[gap.end]-A$Timestamp[gap.start]
  ) 
 ) 
   }
  } 
 }
 x <- B$GapLength
 h<-hist(x, xlab="Gap Length [s?]", 

There's a problem with the rbind function which I'm using wrong.
Is this the right approach and could you please give me a hint on how
to add the line?
Or is there a better way to achieve this?


Best
Ramon


[1]
http://stackoverflow.com/questions/7309411/how-to-calculate-autocorrelation-in-r-zoo-object

[2]
http://www.statmethods.net/graphs/images/histogram3.jpg

__
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] combining two different matrizes

2013-06-06 Thread Berend Hasselman

On 05-06-2013, at 23:56, ThomasH  wrote:

> 
> Hello together,
> 
> this is ma first post, so please aplogize me if post this in the wrong
> section.
> 
> I have problem concerning ma two matrizes.
> 
> After a regressione and so on, I got two matrizes
> 
> Matrixres contains the results of ma calculation.
> 
> Matrixr contains my detiene, which where Aldo used for the regression.
> 
> Please ser the following code:
> 
> #Datei einlesen
> residual = read.csv2("E:***Input-R_Renditen.csv",header=TRUE, sep=";")
> 
> 
> #Aktientitel
> alist <- list()
> for (a in 2:11){
> 
> 
> #Regression
>#Länge Gesamtzeit
>t <- 243
>tx <- t-59
> 
>#Länge Regression
>reglist <- list()
>for (i in 1:tx){
>j <- i+59
> 
>#RegressionsVariable
>x = residual[i:j,a]
>rm = residual[i:j,12]
>smb = residual[i:j,13]
>hml = residual[i:j,14]
>rf = residual[i:j,15]
> 
>#Überschussrenditen
>ex=x-rf
>erm=rm-rf
> 
>#Regression
>reg <- lm(ex~erm+smb+hml)
>reglist[[i]] <- coef(reg)
>
> 
> #Berechnung Residuum
>   #Residual Berechnung
>   rx = residual[(j-5):j,a]
>   rrm = residual[(j-5):j,12]
>   rsmb = residual[(j-5):j,13]
>   rhml = residual[(j-5):j,14]
>   rrf = residual[(j-5):j,15]
> 
>   rex = rx-rrf
>   rerm = rrm-rrf
> 
>   #Berechnung
>   res <-
> sum(rex-(reglist[[i]][2]*rerm+reglist[[i]][3]*rsmb+reglist[[i]][4]*rhml))/sd(rex-(reglist[[i]][2]*rerm+reglist[[i]][3]*rsmb+reglist[[i]][4]*rhml))
>   reglist[[i]] <- res
> }   
>
> 
> #Residuen auf alle Aktien
>alist[[a]] <- reglist
> }
>
> #Matrix mit Residuen
>matrixres <- do.call(cbind,alist)
> 
> #Spaltennamen/Zeilennamen
> s<- names(residual)[2:11]
> colnames(matrixres)<-s
> 
> #RenditeMatrix
> matrixr <- do.call(cbind,residual[60:243,2:11])
> 
> 
> Now I want to combines  the two matrizes in the following way:
> 
> Under every row of matrixres should stand the row of matrixr for excample:
> 
> Matrixres row1
> Matrixr row1
> Matrixres row2
> Matrixr row 2
> 
> Can anybody help me? I was working on this problem the whole day, but have
> no idea.

Something like this (assuming A and B have the same dimensions)

set.seed(11)
A <- matrix(round(rnorm(25),3),nrow=5)
B <- matrix(1:25,nrow=5)

C <- matrix(0,nrow=2*nrow(A),ncol=ncol(A))

crows <- seq.int(from=1,to=2*nrow(A),by=2)
C[crows,] <- A
C[crows+1,] <- B

Some friendly advice: get someone to check your English before sending a mail 
to the list.

Berend

__
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] Loop through variables and estimate effects on several outcomes

2013-06-06 Thread arun
Hi,
Try:
hsb2 <- read.csv("http://www.ats.ucla.edu/stat/data/hsb2.csv";)

varlist<-names(hsb2)[8:10]
fun2<- function(varName){
    res<- sapply(varName,function(x){
              model1<- 
lm(substitute(cbind(female,race,ses)~i,list(i=as.name(x))),data=hsb2)
          sM<- summary(model1)
              sapply(sM,function(x) x$coef[2,1])           
             })
            res 
    }     

 fun2(varlist)
# write math  science
#Response female 0.01350896 -0.001563341 -0.006441112
#Response race   0.02412624  0.022474213  0.033622966
#Response ses    0.01585530  0.021064315  0.020692042

A.K.

>This post has NOT been accepted by the mailing list yet. 
>I want to estimate the effects of an exposure on several outcomes. The example 
>in this link provides how to loop though variables which are 
>explanatory variables.  
>http://www.ats.ucla.edu/stat/r/pages/looping_strings.htm
>The
 example below estimates the effects of several variables on read.  But I
 want to estimate the effect of  "female" , "race"  ,  "ses"  on 
 "write" ,  >"math"    "science"   one at a time using the hsb data set. 
 How can I loop through these outcomes? 
>varlist <- names(hsb2)[8:11] 
>models <- lapply(varlist, function(x) { 
 > lm(substitute(read ~ i, list(i = as.name(x))), data = hsb2) 
>})

__
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] Post hoc power analysis for mixed-effects models

2013-06-06 Thread Kota Hattori
Hi David,

Thank you very much for your reply. I do understand your point. I was requested 
to do
power analyses for null findings by reviewers. I think simulations would be an 
alternative
choice given the validity of post hoc power analysis is questionable. The 
package, pamm
does simulation. But, as I mentioned before, pamm seems not to be able to 
handle categorical
variables. Is there any other resource available?

This email may be confidential and subject to legal privilege, it may
not reflect the views of the University of Canterbury, and it is not
guaranteed to be virus free. If you are not an intended recipient,
please notify the sender immediately and erase all copies of the message
and any attachments.

Please refer to http://www.canterbury.ac.nz/emaildisclaimer for more
information.

[[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] Error invalid graphics state using text()

2013-06-06 Thread Daniel Wagner

   Hi all!

   I'm using ssplot for drawing a map of Austria and colour the nine provinces
   regarding their share of employment. Now I wanted to add the figures in each
   province and failed miserably. Using the locator() and text() function
   caused the error message "invalid graphics state". I try to show you what I
   have done below, maybe you can find a general fault in my codes. I know that
   it's probably not the shortest way to the finishing line, but I'm a starter
   and  happy  about  every finish. ;) I'm really hoping for some helpful
   comments!

library(maptools)
library(fields)
library(raster)
library(rgdal)
library(RColorBrewer)
library(rgeos)
library(pixmap)
library(classInt)
library(sp)
NUTS2<-readShapePoly("xxx.shp")
shareub2012 <- read.csv("xxx.csv", header=TRUE, sep = ";", dec = ",",stringsAsF
actors=F)
mergeshareub2012 <- merge(shareub2012, NUTS2, by.x = "Bundesland", by.y = "NAME
" )
sortmergeshare2012 <- mergeshareub2012[order(mergeshareub2012$ID),]
col_no <- as.factor(as.numeric(cut(sortmergeshare2012$x, c(0, .075, .125, .25, 
.50
levels(col_no) <- c("< 7,5%", "7,5-12,5%", "12,5-25%","> 25%")
sortmergeshare2012$col_no <- col_no
NUTS2$col_no <- col_no
myPalette<-brewer.pal(4,"Blues")
shareplot2012 <- spplot(NUTS2, "col_no", col="black", col.regions=myPalette,
par.settings=list(axis.line=list(col="transparent")))
shareplot2012
plot.new()
locator()
text(0.82421448,0.3897917,"12,41", cex=0.7)
__
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] SnpMatrix super slow to access cells when large

2013-06-06 Thread nickDIL
Dear snpStats users,

I'm working with a large SnpMatrix object (roughly 5000 samples x 200K snps)
and I've noticed using numerical accessors is extremely slow, e.g, see times
below, takes over 1.5 seconds to retrieve a single cell in SnpMatrix format
[1,1], versus 0.0 seconds to access the same datapoint in RAW format. It
also takes no longer (still 1.5s) to access an entire row or column [1,] or
[,1].

Is snpStats::SnpMatrix doing something unnecessary prior to returning the
matrix entry? [NB: 'chopsticks' seems to give the same slow result]

Is there any way around this delay other than copying the entire SnpMatrix
into RAW format? I want to access specific cell ranges many times in an
algorithm i'm writing and this would be excessively slow with access times
of 1.5s.

Code to show this below.

Many thanks,

N.

# generate raw matrix
rawd <- as.raw(sample(0:3,(10^9),replace=T)); dim(rawd) <- c(5000,20)

# copy to a SnpMatrix object
snpd <- new("SnpMatrix",rawd)


# show class details
> class(snpd)
[1] "SnpMatrix"
attr(,"package")
[1] "snpStats"


# access times in SnpMatrix format

> system.time(snpd[1,])
   user  system elapsed 
  0.876   0.681   1.554 
> system.time(snpd[1,1])
   user  system elapsed 
  0.872   0.668   1.538 
> system.time(snpd[,1])
   user  system elapsed 
  0.896   0.644   1.540 



# access times in raw format

> system.time(rawd[1,])
   user  system elapsed 
  0.012   0.004   0.011 
> system.time(rawd[,1])
   user  system elapsed 
  0   0   0 
> system.time(rawd[1,1])
   user  system elapsed 
  0   0   0 





--
View this message in context: 
http://r.789695.n4.nabble.com/SnpMatrix-super-slow-to-access-cells-when-large-tp4668812.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] Error invalid graphics state using text()

2013-06-06 Thread danielw7
I'm using ssplot for drawing a map of Austria and colour the nine provinces
regarding their share of employment. Now I wanted to add the figures in each
province and failed miserably. Using the locator() and text() function
caused the error message "invalid graphics state". I try to show you what I
have done below, maybe you can find a general fault in my codes. I know that
it's probably not the shortest way to the finishing line, but I'm a starter
and happy about every finish. ;) I'm really hoping for some helpful
comments!

library(maptools)
library(fields)
library(raster)
library(rgdal)
library(RColorBrewer)
library(rgeos)
library(pixmap)
library(classInt)
library(sp)
NUTS2<-readShapePoly("xxx.shp")
shareub2012 <- read.csv("xxx.csv", header=TRUE, sep = ";", dec =
",",stringsAsFactors=F)
mergeshareub2012 <- merge(shareub2012, NUTS2, by.x = "Bundesland", by.y =
"NAME" )
sortmergeshare2012 <- mergeshareub2012[order(mergeshareub2012$ID),]
col_no <- as.factor(as.numeric(cut(sortmergeshare2012$x, c(0, .075, .125,
.25, .50
levels(col_no) <- c("< 7,5%", "7,5-12,5%", "12,5-25%","> 25%")
sortmergeshare2012$col_no <- col_no
NUTS2$col_no <- col_no
myPalette<-brewer.pal(4,"Blues")
shareplot2012 <- spplot(NUTS2, "col_no", col="black", col.regions=myPalette,
par.settings=list(axis.line=list(col="transparent")))
shareplot2012
plot.new()
locator()
text(0.82421448,0.3897917,"12,41", cex=0.7)



--
View this message in context: 
http://r.789695.n4.nabble.com/Error-invalid-graphics-state-using-text-tp4668809.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] [R-pkgs] bpcp package for pointwise confidence intervals for a survival distribution

2013-06-06 Thread Fay, Michael (NIH/NIAID) [E]
Hi all,

I just uploaded a new version of the bpcp package.  It calculates confidence 
intervals for a survival distribution for right-censored data using the newly 
developed beta product confidence procedure.  Previously developed methods can 
have substantial error rate inflation for the lower limit, especially at the 
right end of the curves when there are small numbers of events.  The bpcp 
method guarantees coverage when there is no censoring (when it reduces to the 
exact Clopper-Pearson intervals for binomial data at each time) or when there 
is progressive type II censoring.  Simulations show at least nominal  coverage 
for independent censoring.  Also, it is fast so it can be used routinely.

For a complete description of the new method see

Fay, Brittain, and Proschan (2013). Pointwise confidence intervals for a 
survival distribution with small samples or heavy censoring. Biostatistics 
doi:10.1093/biostatistics/kxt016 (available at 
http://www.niaid.nih.gov/about/organization/dcr/brb/staff/Pages/michael.aspx ).


Mike

**
Michael P. Fay, PhD
Mathematical Statistician
National Institute of Allergy and Infectious Diseases
Tel: 301-451-5124   Fax:301-480-0912
(U.S. postal mail address)
6700B Rockledge Drive MSC 7609
Bethesda, MD 20892-7609
(Overnight mail address)
6700-A Rockledge Drive, Room 5133
Bethesda, MD 20817
http://www3.niaid.nih.gov/about/organization/dcr/BRB/staff/michael.htm
**
Disclaimer:
The information in this e-mail and any of its attachment...{{dropped:7}}

__
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] Problem with marginal effects of a multinomial logistic regression

2013-06-06 Thread Leonard Moulin
Hi r-users,

I try to calculate marginal effects of a multinomial logistic regression. To do 
this i use mlogit package and effects() function.

Here is how the procedure works (source : effects() function of mlogit package) 
:
data("Fishing", package = "mlogit")
Fish <- mlogit.data(Fishing, varying = c(2:9), shape = "wide", choice = "mode")
m <- mlogit(mode ~ price | income | catch, data = Fish)
# compute a data.frame containing the mean value of the covariates in
# the sample
z <- with(Fish, data.frame(price = tapply(price, index(m)$alt, mean),
   catch = tapply(catch, index(m)$alt, mean),
   income = mean(income)))
# compute the marginal effects (the second one is an elasticity
effects(m, covariate = "income", data = z)
effects(m, covariate = "price", type = "rr", data = z)
effects(m, covariate = "catch", type = "ar", data = z)
I have no problem with first step (mlogit.data() function). I think my problem 
is on the specification of the multinomial regression.

My regression (for example with three variables) is on the form: Y ~ 0 | X1 + 
X2 + X3. When I try to estimate the marginal effects for model 2 variables no 
problem, however for 3 variables R console returns me the following error: 
"Error in if (rhs% in% c (1, 3)) {: argument is of length zero " (translation 
from error in R console in french).

To understand what is my problem I tried to perform a multinomial regression of 
similar shape on the dataset "Fishing", ie: mode ~ 0 | income + price + catch 
(even if this form has no "economic" sense.) Again the R console returns me the 
same error for 3 variables but manages to estimate these effects for a model 
with two variables ...

This leads me to think that my problem really comes from the specification of 
my multinomial regression ... Do you know how I could find a solution to my 
problem?

Thank you for your help :)

Léonard.
[[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 to interval package

2013-06-06 Thread Fay, Michael (NIH/NIAID) [E]
Hi all,

I just uploaded an update to the interval package that does NPMLE estimates of 
survival distribution and weighted logrank tests for interval censored data.

The update now includes a confidence interval method for the survival that uses 
a modified bootstrap method.

Mike

**
Michael P. Fay, PhD
Mathematical Statistician
National Institute of Allergy and Infectious Diseases
Tel: 301-451-5124   Fax:301-480-0912
(U.S. postal mail address)
6700B Rockledge Drive MSC 7609
Bethesda, MD 20892-7609
(Overnight mail address)
6700-A Rockledge Drive, Room 5133
Bethesda, MD 20817
http://www3.niaid.nih.gov/about/organization/dcr/BRB/staff/michael.htm
**
Disclaimer:
The information in this e-mail and any of its attachment...{{dropped:7}}

__
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] basic sweave question

2013-06-06 Thread Michael Liu
Dear Ligges,

*Well, *
*
*
*1st, you should find the "*Sweave.sty" file, I am sure you can.
2nd, in the Latex configure menu, adding the path containing the "
Sweave.sty"

That will be OK.

Best Wishes
  Guanhao Liu


On Fri, Mar 22, 2013 at 6:56 PM, Uwe Ligges  wrote:

> I typically run
>
> R CMD Sweave
>
> which resolves issues with style files in directories  unknown to your
> LaTeX distribution.
>
> Best,
> Uwe Ligges
>
>
> On 21.03.2013 18:44, lgh0504 wrote:
>
>> The pdflatex is the an excutable programe located in your MiKTex‘ bin
>>
>> folder, I guess you have not add your MikTex bin folder into your Windows
>> PATH environment
>>
>>
>> On Thu, Mar 14, 2013 at 12:55 AM, susieboyce [via R] <
>> ml-node+s789695n4661220h61@n4.**nabble.com>
>> wrote:
>>
>>  I have located my Sweave.sty in my R program files and I've tried copying
>>> and pasting this into many folders in the C/.../MiKTeX/tex/latex
>>> environment, including making a new folder called 'Sweave' and storing in
>>> here and still my .tex file gives me an error when I try to compile it.
>>>
>>> What is this pdflatex? Is it an R command? If so, what package is it in?
>>>
>>> Where does this go? In the original .rnw file, do you type this into R or
>>> somewhere else?:
>>> pdflatex --include-directory="C:\**Program
>>> Files\R\R-2.14.0\share\texmf\**tex\latex" your_tex.tex
>>>
>>>
>>> --
>>>
>>>   If you reply to this email, your message will be added to the
>>> discussion
>>> below:
>>> http://r.789695.n4.nabble.com/**basic-sweave-question-**
>>> tp876734p4661220.html
>>>   To unsubscribe from basic sweave question, click here<
>>> http://r.789695.n4.**nabble.com/template/**NamlServlet.jtp?macro=**
>>> unsubscribe_by_code&node=**876734&code=**bGdoMDUwNEBnbWFpbC5jb218ODc2Nz*
>>> *M0fC0xODczMTIwODc1
>>> >
>>> .
>>> NAML>> NamlServlet.jtp?macro=macro_**viewer&id=instant_html%**
>>> 21nabble%3Aemail.naml&base=**nabble.naml.namespaces.**
>>> BasicNamespace-nabble.view.**web.template.NabbleNamespace-**
>>> nabble.view.web.template.**NodeNamespace&breadcrumbs=**
>>> notify_subscribers%21nabble%**3Aemail.naml-instant_emails%**
>>> 21nabble%3Aemail.naml-send_**instant_email%21nabble%**3Aemail.naml
>>> >
>>>
>>>
>>
>>
>>
>> --
>> View this message in context: http://r.789695.n4.nabble.com/**
>> basic-sweave-question-**tp876734p4662100.html
>> Sent from the R help mailing list archive at Nabble.com.
>> [[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] SPlus script

2013-06-06 Thread Scott Raynaud
Ok.  Now I see that the sshc function is not being called.  Thanks for pointing 
that out.  
I'm not certain about the solution, however.  I tried putting call("sshc") at 
the end of the 
program, but nothing happened.  My memory about all of this is fuzzy.  
Suggestions 
on how to call the function appreciated.

- Original Message -
From: William Dunlap 
To: Scott Raynaud ; "r-help@r-project.org" 

Cc: 
Sent: Wednesday, June 5, 2013 2:17 PM
Subject: RE: [R] SPlus script

Both the R and S+ versions (which seem to differ only in the use of _ for 
assignment
in the S+ version) do nothing but define some functions.  You would not expect 
any
printed output unless you used those functions on some data.  Is there another 
script
that does that?

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 Scott Raynaud
> Sent: Wednesday, June 05, 2013 6:21 AM
> To: r-help@r-project.org
> Subject: [R] SPlus script
> 
> This originally was an SPlus script that I modifeid about a year-and-a-half 
> ago.  It worked
> perfectly then.  Now I can't get any output despite not receiving an error 
> message.  I'm
> providing the SPLUS script as a reference.  I'm running R15.2.2.  Any help 
> appreciated.
> 
> MY
> MODIFICATION***
> **
> ## sshc.ssc: sample size calculation for historical control studies
> ## J. Jack Lee (jj...@mdanderson.org) and Chi-hong Tseng
> ## Department of Biostatistics, Univ. of Texas M.D. Anderson Cancer Center
> ##
> ## 3/1/99
> ## updated 6/7/00: add loess
> ##--
>  Required Input:
> #
> # rc number of response in historical control group
> # nc sample size in historical control
> # d  target improvement = Pe - Pc
> # method 1=method based on the randomized design
> #    2=Makuch & Simon method (Makuch RW, Simon RM. Sample size 
> considerations
> #  for non-randomized comparative studies. J of Chron Dis 1980; 
> 3:175-181.
> #    3=uniform power method
>  optional Input:
> #
> # alpha  size of the test
> # power  desired power of the test
> # tol    convergence criterion for methods 1 & 2 in terms of sample size
> # tol1   convergence criterion for method 3 at any given obs Rc in terms of 
> difference
> #  of expected power from target
> # tol2   overall convergence criterion for method 3 as the max absolute 
> deviation
> #  of expected power from target for all Rc
> # cc range of multiplicative constant applied to the initial values ne
> # l.span smoothing constant for loess
> #
> # Note:  rc is required for methods 1 and 2 but not 3
> #    method 3 return the sample size need for rc=0 to (1-d)*nc
> #
>  Output
> # for methdos 1 & 2: return the sample size needed for the experimental group 
> (1
> number)
> #    for given rc, nc, d, alpha, and power
> # for method 3:  return the profile of sample size needed for given nc, 
> d, alpha, and
> power
> #    vector $ne contains the sample size corresponding to 
> rc=0, 1, 2, ... nc*(1-d)
> #    vector $Ep contains the expected power corresponding to
> #  the true pc = (0, 1, 2, ..., nc*(1-d)) / nc
> #
> #--
> sshc<-function(rc, nc=1092, d=.085779816, method=3, alpha=0.05, power=0.8,
>   tol=0.01, tol1=.0001, tol2=.005, cc=c(.1,2), l.span=.5)
> {
> ### for method 1
> if (method==1) {
>  ne1<-ss.rand(rc,nc,d,alpha=.05,power=.8,tol=.01)
>  return(ne=ne1)
>    }
> ### for method 2
> if (method==2) {
> ne<-nc
> ne1<-nc+50
> while(abs(ne-ne1)>tol & ne1<10){
> ne<-ne1
> pe<-d+rc/nc
> ne1<-nef(rc,nc,pe*ne,ne,alpha,power)
> ## if(is.na(ne1)) print(paste('rc=',rc,',nc=',nc,',pe=',pe,',ne=',ne))
> }
> if (ne1>10) return(NA)
> else return(ne=ne1)
> }
> ### for method 3
> if (method==3) {
> if (tol1 > tol2/10) tol1<-tol2/10
> ncstar<-(1-d)*nc
> pc<-(0:ncstar)/nc
> ne<-rep(NA,ncstar + 1)
> for (i in (0:ncstar))
> { ne[i+1]<-ss.rand(i,nc,d,alpha=.05,power=.8,tol=.01)
> }
> plot(pc,ne,type='l',ylim=c(0,max(ne)*1.5))
> ans<-c.searchd(nc, d, ne, alpha, power, cc, tol1)
> ### check overall absolute deviance
> old.abs.dev<-sum(abs(ans$Ep-power))
> ##bad<-0
> print(round(ans$Ep,4))
> print(round(ans$ne,2))
> lines(pc,ans$ne,lty=1,col=8)
> old.ne<-ans$ne
> ##while(max(abs(ans$Ep-power))>tol2 & bad==0){   unnecessary ##
> while(max(abs(ans$Ep-power))>tol2){
> ans<-c.searchd(nc, d, ans$ne, alpha, power, cc, tol1)
> abs.dev<-sum(abs(ans$Ep-power))
> print(paste(" old.abs.dev=",old.abs.dev))
> print(paste(" abs.dev=",abs.dev))
> ##if (abs.dev > old.abs.dev) { bad<-1}
> old.abs.dev<-abs.dev
> print(round(ans$Ep,4))
> print(round(ans$ne,2))

Re: [R] lme function cannot find object

2013-06-06 Thread Ben Bolker
Pfeiffer, Steven  miamioh.edu> writes:


> I have been using the function lme() from package 'nlme' for several months
> now without any problems.  Suddenly, it cannot find a factor in my data.
> Is this a new bug of some kind?  My code and output are below.
> Thanks for your help!
> -Steve Pfeiffer

  You have an extra comma at the end of the first line of
your lme code (indicated with ^^^), so R thinks your interaction term is a
separate argument and tries to intepret it outside the context
of the formula ...

library("nlme")
fit.1<-lme(soil.moisture ~ Trenching + DaysSinceEvent,
   ^
 + Trenching:DaysSinceEvent,
 random = ~ DaysSinceEvent | Plot,
 data=Dat.middle,
 method="ML"
)

__
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] SPlus script

2013-06-06 Thread Ista Zahn
Hi Scott,

As others have pointed out, all your script does is define functions.
It doesn't generate "output" because those functions are never called.
Basically you are missing part of the program--possibly in another
file--that actually calls the sshc function.

Best,
Ista

On Thu, Jun 6, 2013 at 7:24 AM, Scott Raynaud  wrote:
> Ok. I tried copying the original program, changing all _ to <- and running in 
> 3.0.1.  Still no error
> message or output.  It's a mystery to me.
>
> - Original Message -
> From: William Dunlap 
> To: Scott Raynaud ; "r-help@r-project.org" 
> 
> Cc:
> Sent: Wednesday, June 5, 2013 2:17 PM
> Subject: RE: [R] SPlus script
>
> Both the R and S+ versions (which seem to differ only in the use of _ for 
> assignment
> in the S+ version) do nothing but define some functions.  You would not 
> expect any
> printed output unless you used those functions on some data.  Is there 
> another script
> that does that?
>
> 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 Scott Raynaud
>> Sent: Wednesday, June 05, 2013 6:21 AM
>> To: r-help@r-project.org
>> Subject: [R] SPlus script
>>
>> This originally was an SPlus script that I modifeid about a year-and-a-half 
>> ago.  It worked
>> perfectly then.  Now I can't get any output despite not receiving an error 
>> message.  I'm
>> providing the SPLUS script as a reference.  I'm running R15.2.2.  Any help 
>> appreciated.
>>
>> MY
>> MODIFICATION***
>> **
>> ## sshc.ssc: sample size calculation for historical control studies
>> ## J. Jack Lee (jj...@mdanderson.org) and Chi-hong Tseng
>> ## Department of Biostatistics, Univ. of Texas M.D. Anderson Cancer Center
>> ##
>> ## 3/1/99
>> ## updated 6/7/00: add loess
>> ##--
>>  Required Input:
>> #
>> # rc number of response in historical control group
>> # nc sample size in historical control
>> # d  target improvement = Pe - Pc
>> # method 1=method based on the randomized design
>> #2=Makuch & Simon method (Makuch RW, Simon RM. Sample size 
>> considerations
>> #  for non-randomized comparative studies. J of Chron Dis 1980; 
>> 3:175-181.
>> #3=uniform power method
>>  optional Input:
>> #
>> # alpha  size of the test
>> # power  desired power of the test
>> # tolconvergence criterion for methods 1 & 2 in terms of sample size
>> # tol1   convergence criterion for method 3 at any given obs Rc in terms of 
>> difference
>> #  of expected power from target
>> # tol2   overall convergence criterion for method 3 as the max absolute 
>> deviation
>> #  of expected power from target for all Rc
>> # cc range of multiplicative constant applied to the initial values ne
>> # l.span smoothing constant for loess
>> #
>> # Note:  rc is required for methods 1 and 2 but not 3
>> #method 3 return the sample size need for rc=0 to (1-d)*nc
>> #
>>  Output
>> # for methdos 1 & 2: return the sample size needed for the experimental 
>> group (1
>> number)
>> #for given rc, nc, d, alpha, and power
>> # for method 3:  return the profile of sample size needed for given nc, 
>> d, alpha, and
>> power
>> #vector $ne contains the sample size corresponding to 
>> rc=0, 1, 2, ... nc*(1-d)
>> #vector $Ep contains the expected power corresponding to
>> #  the true pc = (0, 1, 2, ..., nc*(1-d)) / nc
>> #
>> #--
>> sshc<-function(rc, nc=1092, d=.085779816, method=3, alpha=0.05, power=0.8,
>>   tol=0.01, tol1=.0001, tol2=.005, cc=c(.1,2), l.span=.5)
>> {
>> ### for method 1
>> if (method==1) {
>>  ne1<-ss.rand(rc,nc,d,alpha=.05,power=.8,tol=.01)
>>  return(ne=ne1)
>>}
>> ### for method 2
>> if (method==2) {
>> ne<-nc
>> ne1<-nc+50
>> while(abs(ne-ne1)>tol & ne1<10){
>> ne<-ne1
>> pe<-d+rc/nc
>> ne1<-nef(rc,nc,pe*ne,ne,alpha,power)
>> ## if(is.na(ne1)) print(paste('rc=',rc,',nc=',nc,',pe=',pe,',ne=',ne))
>> }
>> if (ne1>10) return(NA)
>> else return(ne=ne1)
>> }
>> ### for method 3
>> if (method==3) {
>> if (tol1 > tol2/10) tol1<-tol2/10
>> ncstar<-(1-d)*nc
>> pc<-(0:ncstar)/nc
>> ne<-rep(NA,ncstar + 1)
>> for (i in (0:ncstar))
>> { ne[i+1]<-ss.rand(i,nc,d,alpha=.05,power=.8,tol=.01)
>> }
>> plot(pc,ne,type='l',ylim=c(0,max(ne)*1.5))
>> ans<-c.searchd(nc, d, ne, alpha, power, cc, tol1)
>> ### check overall absolute deviance
>> old.abs.dev<-sum(abs(ans$Ep-power))
>> ##bad<-0
>> print(round(ans$Ep,4))
>> print(round(ans$ne,2))
>> lines(pc,ans$ne,lty=1,col=8)
>> old.ne<-ans$ne
>> ##while(max(abs(ans$Ep-power))>tol2 & bad==0){   unnecessary ##
>> while(max(ab

Re: [R] SPlus script

2013-06-06 Thread Scott Raynaud
Ok. I tried copying the original program, changing all _ to <- and running in 
3.0.1.  Still no error 
message or output.  It's a mystery to me.

- Original Message -
From: William Dunlap 
To: Scott Raynaud ; "r-help@r-project.org" 

Cc: 
Sent: Wednesday, June 5, 2013 2:17 PM
Subject: RE: [R] SPlus script

Both the R and S+ versions (which seem to differ only in the use of _ for 
assignment
in the S+ version) do nothing but define some functions.  You would not expect 
any
printed output unless you used those functions on some data.  Is there another 
script
that does that?

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 Scott Raynaud
> Sent: Wednesday, June 05, 2013 6:21 AM
> To: r-help@r-project.org
> Subject: [R] SPlus script
> 
> This originally was an SPlus script that I modifeid about a year-and-a-half 
> ago.  It worked
> perfectly then.  Now I can't get any output despite not receiving an error 
> message.  I'm
> providing the SPLUS script as a reference.  I'm running R15.2.2.  Any help 
> appreciated.
> 
> MY
> MODIFICATION***
> **
> ## sshc.ssc: sample size calculation for historical control studies
> ## J. Jack Lee (jj...@mdanderson.org) and Chi-hong Tseng
> ## Department of Biostatistics, Univ. of Texas M.D. Anderson Cancer Center
> ##
> ## 3/1/99
> ## updated 6/7/00: add loess
> ##--
>  Required Input:
> #
> # rc number of response in historical control group
> # nc sample size in historical control
> # d  target improvement = Pe - Pc
> # method 1=method based on the randomized design
> #    2=Makuch & Simon method (Makuch RW, Simon RM. Sample size 
> considerations
> #  for non-randomized comparative studies. J of Chron Dis 1980; 
> 3:175-181.
> #    3=uniform power method
>  optional Input:
> #
> # alpha  size of the test
> # power  desired power of the test
> # tol    convergence criterion for methods 1 & 2 in terms of sample size
> # tol1   convergence criterion for method 3 at any given obs Rc in terms of 
> difference
> #  of expected power from target
> # tol2   overall convergence criterion for method 3 as the max absolute 
> deviation
> #  of expected power from target for all Rc
> # cc range of multiplicative constant applied to the initial values ne
> # l.span smoothing constant for loess
> #
> # Note:  rc is required for methods 1 and 2 but not 3
> #    method 3 return the sample size need for rc=0 to (1-d)*nc
> #
>  Output
> # for methdos 1 & 2: return the sample size needed for the experimental group 
> (1
> number)
> #    for given rc, nc, d, alpha, and power
> # for method 3:  return the profile of sample size needed for given nc, 
> d, alpha, and
> power
> #    vector $ne contains the sample size corresponding to 
> rc=0, 1, 2, ... nc*(1-d)
> #    vector $Ep contains the expected power corresponding to
> #  the true pc = (0, 1, 2, ..., nc*(1-d)) / nc
> #
> #--
> sshc<-function(rc, nc=1092, d=.085779816, method=3, alpha=0.05, power=0.8,
>   tol=0.01, tol1=.0001, tol2=.005, cc=c(.1,2), l.span=.5)
> {
> ### for method 1
> if (method==1) {
>  ne1<-ss.rand(rc,nc,d,alpha=.05,power=.8,tol=.01)
>  return(ne=ne1)
>    }
> ### for method 2
> if (method==2) {
> ne<-nc
> ne1<-nc+50
> while(abs(ne-ne1)>tol & ne1<10){
> ne<-ne1
> pe<-d+rc/nc
> ne1<-nef(rc,nc,pe*ne,ne,alpha,power)
> ## if(is.na(ne1)) print(paste('rc=',rc,',nc=',nc,',pe=',pe,',ne=',ne))
> }
> if (ne1>10) return(NA)
> else return(ne=ne1)
> }
> ### for method 3
> if (method==3) {
> if (tol1 > tol2/10) tol1<-tol2/10
> ncstar<-(1-d)*nc
> pc<-(0:ncstar)/nc
> ne<-rep(NA,ncstar + 1)
> for (i in (0:ncstar))
> { ne[i+1]<-ss.rand(i,nc,d,alpha=.05,power=.8,tol=.01)
> }
> plot(pc,ne,type='l',ylim=c(0,max(ne)*1.5))
> ans<-c.searchd(nc, d, ne, alpha, power, cc, tol1)
> ### check overall absolute deviance
> old.abs.dev<-sum(abs(ans$Ep-power))
> ##bad<-0
> print(round(ans$Ep,4))
> print(round(ans$ne,2))
> lines(pc,ans$ne,lty=1,col=8)
> old.ne<-ans$ne
> ##while(max(abs(ans$Ep-power))>tol2 & bad==0){   unnecessary ##
> while(max(abs(ans$Ep-power))>tol2){
> ans<-c.searchd(nc, d, ans$ne, alpha, power, cc, tol1)
> abs.dev<-sum(abs(ans$Ep-power))
> print(paste(" old.abs.dev=",old.abs.dev))
> print(paste(" abs.dev=",abs.dev))
> ##if (abs.dev > old.abs.dev) { bad<-1}
> old.abs.dev<-abs.dev
> print(round(ans$Ep,4))
> print(round(ans$ne,2))
> lines(pc,old.ne,lty=1,col=1)
> lines(pc,ans$ne,lty=1,col=8)
> ### add convex
> ans$ne<-convex(pc,ans$ne)$wy
> ### add loess
> ###old.ne<-ans$ne
> loess.ne<-loess(ans$

Re: [R] Send Mail R and Socket Connections

2013-06-06 Thread Uwe Ligges



On 05.06.2013 11:56, TwistedSkies wrote:

Good Afternoon All,

I am attempting to use the SendMailR function, I have checked with our I.T.
department that I am using the correct server and I have the right
permissions to connect and they have sent emails via this server but not
through R and I have also checked that the port should be 25.

The code:

/

# E-Mail #

library(sendmailR)

from <- "da...@work.com"
to <- "a...@work.com"
subject <- "Send Mail R- Test"
body <- "TESTING TESTING TESTING"
mailControl=list(smtpServer="uksmtp.global.local")

sendmail(from=from,to=to,subject=subject,msg=body,control=mailControl)
/



I receive the below error:

/
function (host = "localhost", port, server = FALSE, blocking = FALSE,
open = "a+", encoding = getOption("encoding"), timeout =
getOption("timeout"))
.Internal(socketConnection(host, port, server, blocking, open,
encoding, timeout))


/




This is not an error but the definition of the socketConnection() function.

If this is really the result, please contact the package maintainer 
(CCIng) for details, assuming you have the most recent versions of R and 
that package.


Best,
Uwe Ligges





So I figured its an error with or I needed to define a new socket
connection, is this where the problem lies? Could anyone give me any
pointers on where to go next with this to get it working?

Many Thanks in Advance!




--
View this message in context: 
http://r.789695.n4.nabble.com/Send-Mail-R-and-Socket-Connections-tp4668699.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.


Re: [R] multivariate multilevel model

2013-06-06 Thread Jose Iparraguirre
Dear Patty,

I would leave others to point you to the right packages, but I'd suggest you 
read this excellent R-based book: "Data Analysis Using Regression and 
Multilevel/Hierarchical Models", by Andrew Gelman and Jennifer Hill (Cambridge 
University Press, 2006).
Kind regards,

José


Prof. José Iparraguirre
Chief Economist
Age UK

Age UK
Tavis House, 1- 6 Tavistock Square
London, WC1H 9NB

-Original Message-
From: r-help-boun...@r-project.org [mailto:r-help-boun...@r-project.org] On 
Behalf Of Patty Haaem
Sent: 06 June 2013 06:36
To: r-help@r-project.org
Subject: [R] multivariate multilevel model

Hi
I want to fit a multivariate multilevel model for a longitudinal data set. can 
you introduce me a package or some function to fit this model?
Thanks in advance.
[[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.

Age UK Improving later life

www.ageuk.org.uk


 
---
Age UK is a registered charity and company limited by guarantee, (registered 
charity number 1128267, registered company number 6825798). 
Registered office: Tavis House, 1-6 Tavistock Square, London WC1H 9NA.

For the purposes of promoting Age UK Insurance, Age UK is an Appointed 
Representative of Age UK Enterprises Limited, Age UK is an Introducer 
Appointed Representative of JLT Benefit Solutions Limited and Simplyhealth 
Access for the purposes of introducing potential annuity and health 
cash plans customers respectively.  Age UK Enterprises Limited, JLT Benefit 
Solutions Limited and Simplyhealth Access are all authorised and 
regulated by the Financial Services Authority. 
--

This email and any files transmitted with it are confidential and intended 
solely for the use of the individual or entity to whom they are 
addressed. If you receive a message in error, please advise the sender and 
delete immediately.

Except where this email is sent in the usual course of our business, any 
opinions expressed in this email are those of the author and do not 
necessarily reflect the opinions of Age UK or its subsidiaries and associated 
companies. Age UK monitors all e-mail transmissions passing 
through its network and may block or modify mails which are deemed to be 
unsuitable.

Age Concern England (charity number 261794) and Help the Aged (charity number 
272786) and their trading and other associated companies merged 
on 1st April 2009.  Together they have formed the Age UK Group, dedicated to 
improving the lives of people in later life.  The three national 
Age Concerns in Scotland, Northern Ireland and Wales have also merged with Help 
the Aged in these nations to form three registered charities: 
Age Scotland, Age NI, Age Cymru.




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