[R] inverse table

2016-06-15 Thread Patrizio Frederic
Dear R-users,
I've a problem that puzzle me

suppose I have a two way contigency  table

a <- sample(al <- letters[1:10],100,T)
b <- sample(bl <- LETTERS[1:5],100,T)
ab <- cbind(a,b)

ddd <- (xtabs(data = ab))
ddd <- as.matrix(ddd)

the question is: how do I reverse the code, thus how do I get raw data
(object ab) from ddd?

I've tried

as.data.frame.table(ddd)

which is not the answer I'm looking for.
Thanks in advance,

PF



-- 
+---
| Patrizio Frederic,
| http://morgana.unimore.it/frederic_patrizio/
+---

[[alternative HTML version deleted]]

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


Re: [R] inverse table

2016-06-16 Thread Patrizio Frederic
Thank you all,
David's solution is the one that matches my taste the most.

Patrizio

On Wed, Jun 15, 2016 at 8:04 PM, David L Carlson  wrote:

> After converting the table to a data frame, replicate each row by the
> number of observations:
>
> > ddd.df <- as.data.frame(ddd)  # as.data.frame.table does the same thing
> > ddd.new <- as.matrix(ddd.df[rep(seq_along(ddd.df[, 1]), ddd.df$Freq),
> 1:2])
> > head(ddd.new)
> a   b
> 1   "a" "A"
> 1.1 "a" "A"
> 2   "b" "A"
> 2.1 "b" "A"
> 3   "c" "A"
> 4   "d" "A"
> > rownames(ddd.new) <- NULL # Optional - get rid of row names
> > head(ddd.new)
>  a   b
> [1,] "a" "A"
> [2,] "a" "A"
> [3,] "b" "A"
> [4,] "b" "A"
> [5,] "c" "A"
> [6,] "d" "A"
>
> -
> David L Carlson
> Department of Anthropology
> Texas A&M University
> College Station, TX 77840-4352
>
> -Original Message-
> From: R-help [mailto:r-help-boun...@r-project.org] On Behalf Of Patrizio
> Frederic
> Sent: Wednesday, June 15, 2016 11:10 AM
> To: r-help@r-project.org
> Subject: [R] inverse table
>
> Dear R-users,
> I've a problem that puzzle me
>
> suppose I have a two way contigency  table
>
> a <- sample(al <- letters[1:10],100,T)
> b <- sample(bl <- LETTERS[1:5],100,T)
> ab <- cbind(a,b)
>
> ddd <- (xtabs(data = ab))
> ddd <- as.matrix(ddd)
>
> the question is: how do I reverse the code, thus how do I get raw data
> (object ab) from ddd?
>
> I've tried
>
> as.data.frame.table(ddd)
>
> which is not the answer I'm looking for.
> Thanks in advance,
>
> PF
>
>
>
> --
> +---
> | Patrizio Frederic,
> | http://morgana.unimore.it/frederic_patrizio/
> +---
>
> [[alternative HTML version deleted]]
>
> __
> R-help@r-project.org mailing list -- To UNSUBSCRIBE and more, see
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide
> http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
>



-- 
+---
| Patrizio Frederic,
| http://morgana.unimore.it/frederic_patrizio/
+---

[[alternative HTML version deleted]]

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


Re: [R] mle2(): logarithm of negative pdfs

2008-07-23 Thread Patrizio Frederic
this is an off-topic, of course!
by definition, pdf MUST BE POSITIVE and log function is defined ONLY
for positive values.
check your pdf code and find out the errors.

regards

Patrizio Frederic


+-
| Patrizio Frederic
| Research associate in Statistics,
| Department of Economics,
| University of Modena and Reggio Emilia,
| Via Berengario 51,
| 41100 Modena, Italy
|
| tel:  +39 059 205 6727
| fax:  +39 059 205 6947
| mail: [EMAIL PROTECTED]
+-


2008/7/23 Zornitsa Luleva <[EMAIL PROTECTED]>:
> Hi,
>
> In order to use the mle2-function, one has to define the likelihood function
> itself. As we know, the likelihood function is a sum of the logarithm of
> probability density functions (pdf).  I have implemented myself the pdfs
> that I am using. My problem is, that the pdfs values are negative and I
> cann't take the logarithm of them in the log-likelihood function.
>
> So how can one take the logarithm of negative values of the pdfs?
>
> Thanks a lot for your advice!
>
> Zoe
>
>[[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-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 simulate heteroscedasticity (correlation)

2008-07-23 Thread Patrizio Frederic
> Now I also want to generate two correlated variables where the error
> variance vary over the variable-correlation.
> And I want to plot this for showing heteroscedasticity.
>
> Like shown here:
> http://upload.wikimedia.org/wikipedia/de/1/1b/Heteroske2.png
>
> Is that possible with R?
>

of course it is. And it' very simple

seed(123456)
x = rnorm(500,1,1)
b0 = 1 # intercept chosen at your choice
b1 = 1 # coef chosen at your choice
h = function(x) 1+.4*x # h performs heteroscedasticity function (here
I used a linear one)
eps = rnorm(500,0,h(x))
y = b0 + b1*x + eps
plot(x,y)
abline(lsfit(x,y))
abline(b0,b1,col=2)

regards

PF

ps notice that in heteroscedasticity case the random vector (X,Y) is
not a bivariate normal but it is:

Y|X=x ~ normal(b0+b1 x; h(x))

ie every conditional Y is normal

+-
| Patrizio Frederic
| Research associate in Statistics,
| Department of Economics,
| University of Modena and Reggio Emilia,
| Via Berengario 51,
| 41100 Modena, Italy
|
| tel:  +39 059 205 6727
| fax:  +39 059 205 6947
| mail: [EMAIL PROTECTED]
+-

__
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] maximum likelihood method to fit a model

2008-07-23 Thread Patrizio Frederic
dear Silvia,
quoting Venables WN and Ripley DB (1994) Modern Applied Statistics
with S-plus, sringer, pag 185:

"Since explicit expressions for the maximum likelihood estimators are
not usually available estimates MUST be calculate iteratively"

means that glm.fit performs MLE indirectly and efficiently.

Hope it's help

regards

PF


2008/7/23 silvia narduzzi <[EMAIL PROTECTED]>:
> Dear R users,
> I use the glm() function to fit a generalized linear model with gamma 
> distribution function and log link.
> I have read in the help page that the default method used by R is "glm.fit" 
> (iteratively reweighted least squares, IWLS).
> Is it possible to use maximum likelihood method?
>
> Thanks
>
>
>
> Silvia Narduzzi
> Dipartimento di Epidemiologia
> ASL RM E
> Via di S. Costanza, 53
> 00198 Roma
> Tel  +39 06 83060461
> Mail  [EMAIL PROTECTED]
>
> __
> 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] sequential sum of a vector...

2008-07-23 Thread Patrizio Frederic
try this

colSums(matrix(x,8))

regards,

PF


+-
| Patrizio Frederic
| Research associate in Statistics,
| Department of Economics,
| University of Modena and Reggio Emilia,
| Via Berengario 51,
| 41100 Modena, Italy
|
| tel:  +39 059 205 6727
| fax:  +39 059 205 6947
| mail: [EMAIL PROTECTED]
+-


2008/7/23 Shubha Vishwanath Karanth <[EMAIL PROTECTED]>:
> Hi R,
>
>
>
> Let,
>
>
>
> x=1:80
>
>
>
> I want to sum up first 8 elements of x, then again next 8 elements of x,
> then again another 8 elements. So, my new vector should look like:
>
> c(36,100,164,228,292,356,420,484,548,612)
>
>
>
> I used:
>
>
>
> aggregate(x,list(rep(1:10,each=8)),sum)[-1]
>
> or
>
> rowsum(x,group=rep(1:10,each=8))
>
>
>
>
>
> But without grouping, can I achieve the required? Any other ways of
> doing this?
>
>
>
> Thanks, Shubha
>
>
>
> This e-mail may contain confidential and/or privileged i...{{dropped:13}}
>
> __
> 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] Convert list of lists <--> data frame

2008-07-23 Thread Patrizio Frederic
try

?unlist

it may help

regards

+-
| Patrizio Frederic
| Research associate in Statistics,
| Department of Economics,
| University of Modena and Reggio Emilia,
| Via Berengario 51,
| 41100 Modena, Italy
|
| tel:  +39 059 205 6727
| fax:  +39 059 205 6947
| mail: [EMAIL PROTECTED]
+-


2008/7/23 Michael Friendly <[EMAIL PROTECTED]>:
> For a function that takes an argument as a list of lists of parameters, I'd
> like to be able to convert that
> to a data.frame and vice versa, but can't quite figure out how.
>
> pats <- list(structure(list(shape = 0, shape.col = "black", shape.lty = 1,
>   cell.fill = "white", back.fill = "white", label = 1, label.size = 1,
>   ref.col = "gray80", ref.grid = "yes", scale.max = 100), .Names =
> c("shape",
> "shape.col", "shape.lty", "cell.fill", "back.fill", "label",
> "label.size", "ref.col", "ref.grid", "scale.max")), structure(list(
>   shape = 0, shape.col = "black", shape.lty = 1, cell.fill = "pink",
>   back.fill = "white", label = 1, label.size = 1, ref.col = "gray80",
>   ref.grid = "yes", scale.max = 100), .Names = c("shape", "shape.col",
> "shape.lty", "cell.fill", "back.fill", "label", "label.size",
> "ref.col", "ref.grid", "scale.max")), structure(list(shape = 0,
>   shape.col = "black", shape.lty = 1, cell.fill = "red", back.fill =
> "white",
>   label = 1, label.size = 1, ref.col = "gray80", ref.grid = "yes",
>   scale.max = 100), .Names = c("shape", "shape.col", "shape.lty",
> "cell.fill", "back.fill", "label", "label.size", "ref.col", "ref.grid",
> "scale.max")))
>
> So, I want pats.df to have 10 columns,
> c("shape",  "shape.col", "shape.lty", "cell.fill", "back.fill", "label",
> "label.size", "ref.col", "ref.grid", "scale.max"), and 3 rows for this
> example.
>
> Given pats.df, I'd want to turn that back to pats.
>
> thanks for any help,
> -Michael
>
> --
> Michael Friendly Email: friendly AT yorku DOT ca Professor, Psychology
> Dept.
> York University  Voice: 416 736-5115 x66249 Fax: 416 736-5814
> 4700 Keele Streethttp://www.math.yorku.ca/SCS/friendly.html
> Toronto, ONT  M3J 1P3 CANADA
>
> __
> 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] How to delete duplicate cases?

2008-07-24 Thread Patrizio Frederic
this works

cno = c(rep(1342,times=3),rep(2568,times=2))
rank= c(.23,.14,.56,.15,.89)

df1 = data.frame(cno,rank)[order(cno,rank),]

cnou= unique(cno)
ind = match(cno,cnou)
where   = tapply(rank,ind,length)
where   = cumsum(as.numeric(where))

df1[where,]


regards,

PF


2008/7/24 Daniel Wagner <[EMAIL PROTECTED]>:
> Dear R users,
>
> I have a dataframe with lot of duplicate cases and I want to delete duplicate 
> ones which have low rank and keep that case which has highest rank.
> e.g
>
>> df1
>   cno  rank
> 1  13420.23
> 2  13420.14
> 3  13420.56
> 4  25680.15
> 5  25680.89
>
> so I want to keep 3rd and 5th  cases with highest rank (0.56 & 0.89) and 
> delete rest of the duplicate cases.
> Could somebody help me?
>
> Regards
>
> Daniel
> Amsterdam
>
>
>
>
>
>
>
>
>
> Send instant messages to your online friends http://uk.messenger.yahoo.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.
>
>

__
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] product of successive rows

2008-07-28 Thread Patrizio Frederic
this works too:

n   = 6 # number of rows
m   = 4 # number of coloumns
nm  = n*m
mat = matrix(1:nm,n) # your matrix
pf  = function(Col){
ind = rep(1:(n/2),each=2)
out = tapply(Col,ind,prod)
out
}
# pf performs forall vecotr x: x[i]*x[i-1], i=2,4,6,...,n

apply(mat,2,pf)

# apply pf to each coloumn of mat


2008/7/28 jim holtman <[EMAIL PROTECTED]>:
> Does this do what you want:
>
>> x <- matrix(1:36,6)
>> x
> [,1] [,2] [,3] [,4] [,5] [,6]
> [1,]17   13   19   25   31
> [2,]28   14   20   26   32
> [3,]39   15   21   27   33
> [4,]4   10   16   22   28   34
> [5,]5   11   17   23   29   35
> [6,]6   12   18   24   30   36
>> # create indices (going to assume an even number of rows
>> x.ind <- seq(1, nrow(x), by=2)
>> t(sapply(x.ind, function(.ind) x[.ind,] * x[.ind+1,]))
> [,1] [,2] [,3] [,4] [,5] [,6]
> [1,]2   56  182  380  650  992
> [2,]   12   90  240  462  756 1122
> [3,]   30  132  306  552  870 1260
>>
>
>
> On Sun, Jul 27, 2008 at 6:20 PM, rcoder <[EMAIL PROTECTED]> wrote:
>>
>> Hi everyone,
>>
>> I want to perform an operation on a matrx that outputs the product of
>> successive pairs of rows. For example: calculating the product between rows
>> 1 & 2; 3 & 4; 5 & 6...etc.
>>
>> Does anyone know of any readily available functions that can do this?
>>
>> Thanks,
>>
>> rcoder
>>
>>
>> --
>> View this message in context: 
>> http://www.nabble.com/product-of-successive-rows-tp18681259p18681259.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.
>>
>
>
>
> --
> Jim Holtman
> Cincinnati, OH
> +1 513 646 9390
>
> What is the problem you are trying to solve?
>
> __
> 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] tensor product of equi-spaced B-splines in the unit square

2008-07-29 Thread Patrizio Frederic
Dear all,
I need to compute tensor product of B-spline defined over equi-spaced
break-points.
I wrote my own program (it works in a 2-dimensional setting)

library(splines)
# set the break-points
Knots   = seq(-1,1,length=10)
# number of splines
M   = (length(Knots)-4)^2
# short cut to splineDesign function
bspline = function(x) splineDesign(Knots,x,outer.ok = T)
# bivariate tensor product of bspline
btens   = function(x) t(bspline(x[1]))%*%bspline(x[2])
# numebr of points to plot
ng  = 51
# create vectors for plotting
xgr = seq(-1,1,length=ng)
xgr2= expand.grid(xgr,xgr)
# generate random coef. of linear combination
bet   = rnorm(M)
# create matrix for contour-type plot
Bx  = apply(xgr2,1,btens)
Bmat= matrix(t(Bx)%*%bet,ng)
# plot the result
contour(xgr,xgr,Bmat)
persp(xgr,xgr,Bmat,theta=15)

any of you have a better idea (ie more efficient)?
Thanks in advance,

Patrizio Frederic

+-
| Patrizio Frederic
| Research associate in Statistics,
| Department of Economics,
| University of Modena and Reggio Emilia,
| Via Berengario 51,
| 41100 Modena, Italy
|
| tel:  +39 059 205 6727
| fax:  +39 059 205 6947
| mail: [EMAIL PROTECTED]
+-

__
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] table questions

2008-07-29 Thread Patrizio Frederic
# is that what you want?

table(cut(xy,seq(0,max(xy)+.4,by=.4)))

# or this

table(cut(xy,hist(xy)$breaks)) # not the same


regards,

PF

+-
| Patrizio Frederic
| Research associate in Statistics,
| Department of Economics,
| University of Modena and Reggio Emilia,
| Via Berengario 51,
| 41100 Modena, Italy
|
| tel:  +39 059 205 6727
| fax:  +39 059 205 6947
| mail: [EMAIL PROTECTED]
+-



2008/7/29 Edna Bell <[EMAIL PROTECTED]>:
> Hi again!
>
>
> Suppose I have the following:
>
>> xy <- round(rexp(20),1)
>> xy
>  [1] 0.1 3.4 1.6 0.4 1.0 1.4 0.2 0.3 1.6 0.2 0.0 0.1 0.1 1.0 2.0 0.9
> 2.5 0.1 1.5 0.4
>> table(xy)
> xy
>  0 0.1 0.2 0.3 0.4 0.9   1 1.4 1.5 1.6   2 2.5 3.4
>  1   4   2   1   2   1   2   1   1   2   1   1   1
>>
> Is there a way to set things up to have
> 0 - 0.4   0.5 - 0.9  etc. please?
>
> I know there is the cut functions, but breaks are required.  If you
> don't have breaks, what should you do, please?
>
> Would using the breaks from the hist function work appropriately, please?
>
> thanks
> Edna Bell
>
> __
> 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] What's the best way to operate R from within Excel?

2008-07-30 Thread Patrizio Frederic
have a look at

http://sunsite.univie.ac.at/rcom/server/doc/RExcel.html

regards,

PF

2008/7/30 losemind <[EMAIL PROTECTED]>:
>
> Hi all,
>
> How do I operate R from within Excel?
>
> I mean, highlight a bunch of cells, and send to R, and do some statistics in
> R, and return back the numbers and do some plots in R.
>
> For example, I have some parameters for Gaussian Random variable, and I
> would find the most convenient way to send these numbers to R and do some
> plots...
>
> How to do that?
>
> Thanks a lot!
>
> --
> View this message in context: 
> http://www.nabble.com/What%27s-the-best-way-to-operate-R-from-within-Excel--tp18722737p18722737.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] Howto Draw Bimodal Gamma Curve with User Supplied Parameters

2008-07-30 Thread Patrizio Frederic
can you write the bimodal gamma pdf?
if so create your own pdf:

dbigamma=function(x,alpha1,alpha2,beta1,beta2)
{
# ... the bimodal gamma pdf (I can't find it)
}

hist(x,probability=T)
curve(dbigamma(x,alpha1,alpha2,beta1,beta2))

and look at the result

regards,

PF

2008/7/29 Gundala Viswanath <[EMAIL PROTECTED]>:
> Hi,
>
> Suppose I have the following vector (data points):
>> x
>  [1]  36.0  57.3  73.3  92.0 300.4  80.9  19.8  31.4  85.8  44.9  24.6  48.0
>  [13]  28.0  38.3  85.2 103.6 154.4 128.5  38.3  72.4 122.7 123.1  41.8  21.7
>  [25] 143.6 120.2  46.6  29.2  44.8  25.0  57.3  96.4  29.4  62.9  66.4  30.0
>  [37]  24.1  14.8  56.6 102.4 117.5  90.4  37.2  79.6  27.8  17.1  26.6  16.3
>  [49]  41.4  48.9  24.1  23.3   9.9  11.5  15.0  23.6  29.3  27.0  19.2  18.7
>  [61]   4.1  13.0   3.3   5.5  38.2   8.5  39.6  39.2  16.1  35.3  23.3  31.5
>  [73]  38.8  51.5  28.4  18.8  24.1  25.4  28.8  32.8  31.0  28.8  33.3  55.5
>  [85]  39.2  21.0  43.7  16.3  50.6  34.6  66.3  50.5  59.4  46.7  51.9 125.6
>  [97]  69.8  43.7  86.8  50.6 132.4  56.0   6.1   4.9   7.1   7.1  12.8  12.1
> [109] 164.2  69.3  15.6  11.4  34.3   9.2  17.6  21.7  19.2  30.7  61.1  35.8
> [121] 185.8 118.4  13.0   9.6  19.1  45.2  94.5 248.0  56.3  24.4  13.8  12.8
> [133]  35.0  31.6  22.5  50.1  18.7  22.1  28.3  39.5  48.2  33.1  43.5  35.1
> [145]  37.4  30.3  15.8  13.9  15.3  16.1  12.7  11.4  13.0  13.8  31.5  25.3
> [157]  65.2  39.5
>
>
> And the following parameter set (2 component) for gamma function.
>
> comp.1comp.2
> alpha (shape)  2.855444  2.152056
> beta  (scale)   10.418785 39.296224
>
> These params are predefined/precalculated by user.
>
> My question is how can I create a bimodal gamma curve - based on the
> two parameter set -
> on top of the histogram of data points above?
>
> - Gundala Viswanath
> Jakarta - Indonesia
>
> __
> 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] tensor product of equi-spaced B-splines in the unit square

2008-07-31 Thread Patrizio Frederic
dear all,
I apologize for a second post on the same subject. I still have the
problem. I'm going to describe the problem in a different setting:
I have 3 matrix A,B,W such that

dim(A) = c(n,M)
dim(B) = c(n,M)
dim(W) = c(M,M)

what I'm searching for is an efficient computation of vector R,
length(R)=n, where

R[i] = \sum_{j=1}^M \sum_{k=1}^M A[i,j]*B[i,k]*W[j,k] # I apologize
for a bit of LaTex code.

I understand this is about tensor products but I can't figure how to
use neither package tensor nor tensorA (I was not trained in tensor
algebra and right now it's too hot in Italy to me for studing a new
topic).

Any suggestion is welcome.

Regards,

Patrizio Frederic

 +-----
 | Patrizio Frederic
 | Research associate in Statistics,
 | Department of Economics,
 | University of Modena and Reggio Emilia,
 | Via Berengario 51,
 | 41100 Modena, Italy
 |
 | tel:  +39 059 205 6727
 | fax:  +39 059 205 6947
 | mail: [EMAIL PROTECTED]
 +-----


2008/7/29 Patrizio Frederic <[EMAIL PROTECTED]>:
> Dear all,
> I need to compute tensor product of B-spline defined over equi-spaced
> break-points.
> I wrote my own program (it works in a 2-dimensional setting)
>
> library(splines)
> # set the break-points
> Knots   = seq(-1,1,length=10)
> # number of splines
> M   = (length(Knots)-4)^2
> # short cut to splineDesign function
> bspline = function(x) splineDesign(Knots,x,outer.ok = T)
> # bivariate tensor product of bspline
> btens   = function(x) t(bspline(x[1]))%*%bspline(x[2])
> # numebr of points to plot
> ng  = 51
> # create vectors for plotting
> xgr = seq(-1,1,length=ng)
> xgr2= expand.grid(xgr,xgr)
> # generate random coef. of linear combination
> bet   = rnorm(M)
> # create matrix for contour-type plot
> Bx  = apply(xgr2,1,btens)
> Bmat= matrix(t(Bx)%*%bet,ng)
> # plot the result
> contour(xgr,xgr,Bmat)
> persp(xgr,xgr,Bmat,theta=15)
>
> any of you have a better idea (ie more efficient)?
> Thanks in advance,
>
> Patrizio Frederic
>
> +-
> | Patrizio Frederic
> | Research associate in Statistics,
> | Department of Economics,
> | University of Modena and Reggio Emilia,
> | Via Berengario 51,
> | 41100 Modena, Italy
> |
> | tel:  +39 059 205 6727
> | fax:  +39 059 205 6947
> | mail: [EMAIL PROTECTED]
> +-
>

__
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] Random number generation

2008-07-31 Thread Patrizio Frederic
?write.table

could help

PF

2008/7/31 Max <[EMAIL PROTECTED]>:
> Marc,
>
> this is very handy. My next question is, do you know a quick and easy way to
> transfer all of the output to a txt file? (or .xls)?
>
> Thanks,
>
> -Max
>
>
> Marc Schwartz explained on 07/31/2008 :
>>
>> on 07/31/2008 12:24 PM Max wrote:
>>>
>>> Hi Everyone,
>>>
>>> I did a quick search of the list and it looks like this may not have been
>>> asked before... I'm trying to generate a matrix of random numbers between 0
>>> and 1, with 6 columns, 1 rows. About all I know is that runif(1) gives
>>> me the random number I'm looking for.
>>>
>>> Any help would be great!
>>>
>>> thanks,
>>
>> MAT <- matrix(runif(1 * 6), 1, 6)
>>
>>  > str(MAT)
>>  num [1:1, 1:6] 0.753 0.600 0.806 0.713 0.796 ...
>>
>>  > head(MAT, 10)
>> [,1]  [,2]   [,3]   [,4]   [,5]   [,6]
>>  [1,] 0.75343430 0.4993896 0.68554749 0.01924549 0.90579982 0.99606191
>>  [2,] 0.59957219 0.4075650 0.57851744 0.97208426 0.32137505 0.02089689
>>  [3,] 0.80567935 0.5746030 0.16520072 0.92615138 0.01628994 0.90075333
>>  [4,] 0.71270574 0.3252210 0.53765089 0.58930899 0.03053356 0.23282879
>>  [5,] 0.79603691 0.5591622 0.97308348 0.52744458 0.76403708 0.22268021
>>  [6,] 0.49624259 0.5106604 0.06687444 0.48659150 0.29803454 0.91760758
>>  [7,] 0.32921909 0.7784539 0.20468873 0.86730697 0.42581735 0.59344279
>>  [8,] 0.93646405 0.4819996 0.79033546 0.68441917 0.28566573 0.97244395
>>  [9,] 0.02964297 0.5489500 0.64355067 0.87131530 0.58505804 0.06972828
>> [10,] 0.55956266 0.8376349 0.11850374 0.37687892 0.71220844 0.97784727
>>
>>
>>
>> The first argument to runif() is how many random deviates you want to
>> generate. If you need to be able to reproduce the exact sequence again in
>> the future, see ?set.seed.
>>
>> HTH,
>>
>> 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.
>
> __
> 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] contour lines in windows device but neither in pdf nor in postscript

2008-08-01 Thread Patrizio Frederic
library(mvtnorm)
x = seq(-4,4,length=201)
xy= expand.grid(x,x)
sigma = (diag(c(1,1))+1)/2
d2= matrix(dmvnorm(xy,sigma=sigma),201)
xsamp = rmvnorm(200,sigma=sigma)

contour(x,x,d2)
points(xsamp,col=3,pch=16)

pdf("pdftry.pdf")
contour(x,x,d2)
points(xsamp,col=3,pch=16)
dev.off()

postscript("pstry.ps")
contour(x,x,d2)
points(xsamp,col=3,pch=16)
dev.off()

# I can see contour lines in a window device but I can't see them in
files pdftry.pdf and pstry.ps
> version
   _
platform   i386-pc-mingw32
arch   i386
os mingw32
system i386, mingw32
status
major  2
minor  7.1
year   2008
month  06
day23
svn rev45970
language   R
version.string R version 2.7.1 (2008-06-23)

what's going wrong?

Thanks in advance,
Regards,

Patrizio Frederic

__
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] contour lines in windows device but neither in pdf nor in postscript

2008-08-01 Thread Patrizio Frederic
dear Ripley and Diffort,
thank you for the quick reply.
I figured out it was my mistake: I wrote this to the list:

contour(x,x,d2)

in fact my framework I used

contour(x,x,d2,labex=0) # that produced the error

then I learned by myself that if I want to suppress labels I have to use

contour(x,x,d2,drawlabel=F)

and everything works fine now. Thank you again,

Patrizio



2008/8/1 Prof Brian Ripley <[EMAIL PROTECTED]>:
> What viewers are you using?
>
> Works for me (using ghostscript 4.61 and acroread 8.1.2) -- I even tried
> 2.7.1 (as well as R-patched).
>
> On Fri, 1 Aug 2008, Patrizio Frederic wrote:
>
>> library(mvtnorm)
>> x = seq(-4,4,length=201)
>> xy= expand.grid(x,x)
>> sigma = (diag(c(1,1))+1)/2
>> d2= matrix(dmvnorm(xy,sigma=sigma),201)
>> xsamp = rmvnorm(200,sigma=sigma)
>>
>> contour(x,x,d2)
>> points(xsamp,col=3,pch=16)
>>
>> pdf("pdftry.pdf")
>> contour(x,x,d2)
>> points(xsamp,col=3,pch=16)
>> dev.off()
>>
>> postscript("pstry.ps")
>> contour(x,x,d2)
>> points(xsamp,col=3,pch=16)
>> dev.off()
>>
>> # I can see contour lines in a window device but I can't see them in
>> files pdftry.pdf and pstry.ps
>>>
>>> version
>>
>>  _
>> platform   i386-pc-mingw32
>> arch   i386
>> os mingw32
>> system i386, mingw32
>> status
>> major  2
>> minor  7.1
>> year   2008
>> month  06
>> day23
>> svn rev45970
>> language   R
>> version.string R version 2.7.1 (2008-06-23)
>>
>> what's going wrong?
>>
>> Thanks in advance,
>> Regards,
>>
>> Patrizio Frederic
>>
>> __
>> 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.
>>
>
> --
> Brian D. Ripley,  [EMAIL PROTECTED]
> Professor of Applied Statistics,  http://www.stats.ox.ac.uk/~ripley/
> University of Oxford, Tel:  +44 1865 272861 (self)
> 1 South Parks Road, +44 1865 272866 (PA)
> Oxford OX1 3TG, UKFax:  +44 1865 272595
>

__
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 regression with constraints

2008-08-08 Thread Patrizio Frederic
Hi Zhang ,
take a look to sur package

http://www.systemfit.org/

regards,

Patrizio Frederic

+-
| Patrizio Frederic
| Research associate in Statistics,
| Department of Economics,
| University of Modena and Reggio Emilia,
| Via Berengario 51,
| 41100 Modena, Italy
|
| tel:  +39 059 205 6727
| fax:  +39 059 205 6947
| mail: [EMAIL PROTECTED]
+-


2008/8/8 Zhang Yanwei - Princeton-MRAm <[EMAIL PROTECTED]>:
> Hi all,
>   I am running a bivariate regression with the following:
>
> p1=c(184,155,676,67,922,22,76,24,39)
> p2=c(1845,1483,2287,367,1693,488,435,1782,745)
> I1=c(1530,1505,2505,204,2285,269,1271,298,2023)
> I2=c(8238,6247,6150,2748,4361,5549,2657,3533,5415)
> R1=I1-p1
> R2=I2-p2
>
> x1=cbind(p1,R1)
> y1=cbind(p2,R2)
>
> fit1=lm(y1~-1+x1)
> summary(fit1)
>
> Response 2:
> Coefficients:
> Estimate Std. Error t value Pr(>|t|)
> x1p1  -1.4969 2.7004  -0.554  0.59662
> x1R1   3.0937 0.8366   3.698  0.00767 **
>
>
> One can see that in the second regression, i.e. R2~-1+p1+R1,  the coefficient 
> for p1 is not significant. I wonder if I can run this bivariate regression 
> again with the constraint that the coefficient for p1 in the second 
> regression equation  is zero? Thanks a lot.
>
> Sincerely,
> Yanwei Zhang
> Department of Actuarial Research and Modeling
> Munich Re America
> Tel: 609-275-2176
> Email: [EMAIL PROTECTED]<mailto:[EMAIL PROTECTED]>
>
>
>[[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-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] Replace NAs in a range of data frame columns

2009-12-12 Thread Patrizio Frederic
hi Michael,
the following code should work

b <- a[match('first',names(a)): match('last',names(a))]
b[is.na(b)]<-0
a[match('first',names(a)): match('last',names(a))] <- b

cheers,

Patrizio

2009/12/13 Michael Scharkow :
> Dear all,
>
> I'm stuck in a seemingly trivial task that I need to perform for many
> datasets. Basically, I want to replace NA with 0 in a specified range of
> columns in a dataframe. I know the first and last column to be recoded only
> by its name.
>
> I can select the columns starting like this
> a[match('first',names(a)): match('last',names(a))]
>
> The question is how can replace all NA with 0 in this subset of the data?
>
> Thanks and greetings,
> Michael
>
> __
> 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.
>



-- 
+-
| Patrizio Frederic, PhD
| Assistant Professor,
| Department of Economics,
| University of Modena and Reggio Emilia,
| Via Berengario 51,
| 41100 Modena, Italy
|
| tel:  +39 059 205 6727
| fax:  +39 059 205 6947
| mail: patrizio.frede...@unimore.it
+-

__
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] data transformation using gamma

2009-05-07 Thread Patrizio Frederic
Roslina,
this code performs what you need:

dt  = matrix((1:(58*12))/58/12,58) # some numbers
# if dt is a data.frame use dt = as.matrix(dt)
a   = (1:12)/12 # some a coef
b   = (12:1)/12 # some b coef
dtgam   = matrix(pgamma(dt,a,b),58)
# dtgam is the transformation you're looking for

no loop needed no transform function involved
cheers,

Patrizio


2009/5/7 Roslina Zakaria :
> Hi R-users,
>
> I have this code to uniformise the data using gamma:
>
>> length(dp1)
> [1] 696
>> dim(dp1)
> [1] 58 12
>> dim(ahall)
> [1]  1 12
>> dim(bhall)
> [1]  1 12
>
>> trans_dt <- function(dt,a,b)
> + { n1 <- ncol(dt)
> +   n2 <- length(dt)
> +   trans  <- vector(mode='numeric', length=n2)
> +   dim(trans) <- dim(dt)
> +   for (i in 1:n1)
> +   {  dt[,i] <- as.vector(dt[,i])
> +  trans[,i] <- transform(dti,newdt=pgamma(dti,shape= 
> a[1,i],scale=b[1,i])) }
> +   trans
> + }
>
>> trans_dt(dp1,ahall,bhall)
> Error in transform(dti, newdt = pgamma(dti, shape = a[1, i], scale = b[1,  :
>   object "dti" not found
>
> and also try
> trans_dt <- function(dt,a,b)
> { n1 <- ncol(dt)
>   n2 <- length(dt)
>   trans  <- vector(mode='numeric', length=n2)
>   dim(trans) <- dim(dt)
>   for (i in 1:n1)
>   {  dti <- dt[,i]
>  ai  <- a[1,i]
>  bi  <- b[1,i]
>  trans[,i] <- transform(dti,newdt=pgamma(dti,shape= ai,scale=bi)) }
>   trans
> }
>
> trans_dt(dp1,ahall,bhall)
> Error in pgamma(dti, shape = ai, scale = bi) : object "dti" not found
>
>
> Thank you for any help given.
>
>
>
>        [[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.
>
>



-- 
+-
| Patrizio Frederic, PhD
| Assistant Professor,
| Department of Economics,
| University of Modena and Reggio Emilia,
| Via Berengario 51,
| 41100 Modena, Italy
|
| tel:  +39 059 205 6727
| fax:  +39 059 205 6947
| mail: patrizio.frede...@unimore.it
+-

__
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] ReadItem: unknown type 136, perhaps written by later version of R

2009-06-08 Thread Patrizio Frederic
Dear all,
I use at least two pc to perform my data analysis. Both are powered
with R updated at the latest release (currently 2.9.0 under windows xp
pro). I bring .Rdata on my portable drive and use them on any of my pc
(this worked also under (k)ubuntu equipped machines). Please notice
that not all my pc have the same libraries installed. If a not
installed library was called in .First() I got an warning message but
the remaining objects were loaded anyhow.

All that worked straight up to 2.8 series. Some problems arises with
2.9.0 version:

1).RData created with pc1 were not loaded in pc2 (pc2 cant'find
appropriate library to read .RData, I fixed that problem by installing
the appropriate package in pc2). Can't remember right now the explicit
error message, but I remember the specific package was mentioned.

2).RData created with pc2 is not loaded in pc1. I got this message:

Error in load("K:\\.RData") :
  ReadItem: unknown type 136, perhaps written by later version of R

I've no idea of what type 136 is. Any hints?
Thanks in advance,

Patrizio

-- 
+-----
| Patrizio Frederic, PhD
| Assistant Professor,
| Department of Economics,
| University of Modena and Reggio Emilia,
| Via Berengario 51,
| 41100 Modena, Italy
|
| tel:  +39 059 205 6727
| fax:  +39 059 205 6947
| mail: patrizio.frede...@unimore.it
+-

__
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] Interpolate x from y

2009-03-25 Thread Patrizio Frederic
Greg,
it seems an obvious behavior to me

y=c(2,2,2,3,3,3,1)
x=1:length(y)
plot(x,y)
lines(x,approxfun(x,y)(x)) # for every x it exists one only value of y

plot(y,x)
lines(sort(y),approxfun(y,x)(sort(y))) # for some y it exists more
than one value of x!

approxfun return a function. By definition a function maps one value
of domain into one only value of codomain. Would you expect that one y
value returns more than one x? I don't
Be more specific and maybe we can help you

Patrizio

2009/3/25 Greg :
> Is it possible to interpolate a value for x with knowledge of y?
>
> For example, approx(x, y, xout) will give me y's given a set of x's,
> which is opposite to what I'm after.  I've tried switching x and y,
> e.g., approx(y, x, xout), but in a real data set it is possible to
> have more than one y for a given x causing approx() to remove
> coordinates.
>
> Thanks for your help,
>
> Greg.
>
> __
> 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] Converting a Matrix to a Vector

2009-03-25 Thread Patrizio Frederic
try also

m <- matrix( runif(5^2), nrow=5, dimnames = Names<- list(
c("A","B","C","D","E"),
 c("O","P","Q","R","S") ) )

data.frame(expand.grid(Names[[1]],Names[[2]]),as.numeric(m))
data.frame(code=outer(Names[[1]],Names[[2]],paste,sep=".")[1:25],num=as.numeric(m))


Patrizio


2009/3/25 jim holtman :
> Use the 'reshape' package:
>
>> library(reshape)
>> melt(m)
>   X1 X2      value
> 1   A  O 0.26550866
> 2   B  O 0.37212390
> 3   C  O 0.57285336
> 4   D  O 0.90820779
> 5   E  O 0.20168193
> 6   A  P 0.89838968
> 7   B  P 0.94467527
> 8   C  P 0.66079779
> 9   D  P 0.62911404
> 10  E  P 0.06178627
> 11  A  Q 0.20597457
> 12  B  Q 0.17655675
> 13  C  Q 0.68702285
> 14  D  Q 0.38410372
> 15  E  Q 0.76984142
> 16  A  R 0.49769924
> 17  B  R 0.71761851
> 18  C  R 0.99190609
> 19  D  R 0.38003518
> 20  E  R 0.77744522
> 21  A  S 0.93470523
> 22  B  S 0.21214252
> 23  C  S 0.65167377
> 24  D  S 0.1210
> 25  E  S 0.26722067
>>
>
>
> On Wed, Mar 25, 2009 at 3:13 AM, Ken-JP  wrote:
>>
>> Say I have:
>>
>>> set.seed( 1 )
>>> m <- matrix( runif(5^2), nrow=5, dimnames = list( c("A","B","C","D","E"),
>>> c("O","P","Q","R","S") ) )
>>> m
>>          O          P         Q         R         S
>> A 0.2655087 0.89838968 0.2059746 0.4976992 0.9347052
>> B 0.3721239 0.94467527 0.1765568 0.7176185 0.2121425
>> C 0.5728534 0.66079779 0.6870228 0.9919061 0.6516738
>> D 0.9082078 0.62911404 0.3841037 0.3800352 0.121
>> E 0.2016819 0.06178627 0.7698414 0.7774452 0.2672207
>>
>> ---
>>
>> I want to create a vector v from matrix m that looks like this:
>>
>> A.O 0.2655087
>> B.O 0.3721239
>>
>> v <- as.vector( m ) almost gives me what I want, but then I need to take
>> combinations of colnames( m ) and rownames( m ) to get my labels and hope
>> they match up in order: if not, manipulate the order.  This approach feels
>> kludgy...
>>
>> Is this the right approach or is there a better way?
>>
>>
>>
>> --
>> View this message in context: 
>> http://www.nabble.com/Converting-a-Matrix-to-a-Vector-tp22696267p22696267.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.
>>
>
>
>
> --
> Jim Holtman
> Cincinnati, OH
> +1 513 646 9390
>
> What is the problem that you are trying to solve?
>
> __
> 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] same value in column-->delete

2009-03-26 Thread Patrizio Frederic
this works

which.is.not.unique <- apply(x,2,function(x)ifelse(length(unique(x))==1,F,T))
x[,which.is.not.unique]

patrizio

2009/3/26 Duijvesteijn, Naomi :
>
>   Hi Readers,
>
>
>   I have a question.
>
>
>   I have a large dataset and want to throw away columns that have the same
>   value in the column itself and I want to know which column this was.
>
>
>   For example
>
>   > x<-data.frame(id=c(1,2,3), snp1=c("A","G",
>   "G"),snp2=c("G","G","G"),snp3=c("G","G","A"))
>
>   > x
>
>     id snp1 snp2 snp3
>
>   1  1    A    G    G
>
>   2  2    G    G    G
>
>   3  3    G    G    A
>
>
>   Now I want to know that snp2 in monomorphic (the same value for the column)
>   and after I know which column it is I want to take these columns out.
>
>
>   Thanks,
>
>   Naomi
>
>
>
>
>
>   Disclaimer:  De  informatie opgenomen in dit bericht (en bijlagen) kan
>   vertrouwelijk zijn en is uitsluitend bestemd voor de geadresseerde(n).
>   Indien u dit bericht ten onrechte ontvangt, wordt u geacht de inhoud niet te
>   gebruiken, de afzender direct te informeren en het bericht te vernietigen.
>   Aan dit bericht kunnen geen rechten of plichten worden ontleend.
>
>   
>   
>
>   Disclaimer: The information contained in this message may be confidential
>   and is intended to be exclusively for the addressee. Should you receive this
>   message unintentionally, you are expected not to use the contents herein, to
>   notify the sender immediately and to destroy the message. No rights can be
>   derived from this message.
>
> __
> 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] how to get all iterations if I meet NaN?

2009-03-27 Thread Patrizio Frederic
2009/3/27 huiming song :
> hi, everybody, please help me with this question:
>
> If I want to do iteration for 1000 times, however, for the 500th iteration,
> there is NaN appears. Then the iteration will stop. If I don't want the stop
> and want the all the 1000 iterations be done. What shall I do?
>
>
> suppose I have x[1:1000] and z[1:1000],I want to do some calculation for all
> x[1] to x[1000].
>
> z=rep(0,1000)
> for (i in 1:1000){
>  z[i]=sin(1/x[i])
> }
>
> if x[900] is 0, in the above code it will not stop when NaN appears. Suppose
> when sin(1/x[900]) is NaN appears and the iteration will now fulfill the
> rest 100 iterations. How can I write a code to let all the 1000 iterations
> be done?
>

not sure I properly understood. Consider:

 x = seq(-pi,pi,length=1001)
 z = sin(1/x)
# Warning message:
# In sin(1/x) : NaNs produced
x[500:502]; z[500:502]
# [1] -0.006283185  0.0  0.006283185
# [1] -0.8754095NaN  0.8754095

one NaN and one warning have been created, the remaining 1000
calculations has been executed.

lim(x->0)sin(1/x) not exists so sin(1/0) is not a number nan

z1 = z[!is.nan(z)]
x1 = x[!is.nan(z)]

# x and z without the z's nan position

hope that help

Patrizio

__
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] Selecting all rows of factors which have at least one positive value?

2009-04-02 Thread Patrizio Frederic
or the exactly equivalent form:

x[x$X1 %in% unique(x[x$X2>0,"X1"]), ]

Patrizio

2009/4/2 Nutter, Benjamin :
> x <-
> data.frame(matrix(c(rep(11,4),rep(12,3),rep(13,3),rep(0,3),1,rep(0,4),re
> p(1,2)),ncol=2))
>
> id.keep <- unique(subset(x,X2>0)$X1)
>
> x2 <- subset(x,X1 %in% id.keep)
>
> x2
>
> -Original Message-
> From: r-help-boun...@r-project.org [mailto:r-help-boun...@r-project.org]
> On Behalf Of Stephan Lindner
> Sent: Thursday, April 02, 2009 11:26 AM
> To: r-h...@stat.math.ethz.ch
> Subject: [R] Selecting all rows of factors which have at least one
> positive value?
>
> Dear all,
>
> I'm trying to select from a dataframe all rows which correspond to a
> factor (the id variable) for which there exists at least one positive
> value of a certain variable. As an example:
>
> x <-
> data.frame(matrix(c(rep(11,4),rep(12,3),rep(13,3),rep(0,3),1,rep(0,4),re
> p(1,2)),ncol=2))
>
>> x
>
>
>   X1 X2
> 1  11  0
> 2  11  0
> 3  11  0
> 4  11  1
> 5  12  0
> 6  12  0
> 7  12  0
> 8  13  0
> 9  13  1
> 10 13  1
>
>
> and I want to select all rows pertaining to factor levels of X1 for
> which exists at least one "1" for X2. To be clear, I want rows 1:4
> (since there exists at least one observation for X1==11 for which
> X2==1) and rows 8:10 (likewise).
>
> It is easy to obtain the corresponding factor levels (i.e.,
> unique(x$X1[x$X2==1])), but I got stalled selecting the corresponding
> rows. I tried grep, but then I have to loop and concatenate the
> resulting vector. Any ideas?
>
>
> Thanks a lot!
>
>
>        Stephan
>
>
>
>
>
> --
> ---
> Stephan Lindner
> University of Michigan
>
> __
> 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.
>
>
> ===
>
> P Please consider the environment before printing this e-mail
>
> Cleveland Clinic is ranked one of the top hospitals
> in America by U.S. News & World Report (2008).
> Visit us online at http://www.clevelandclinic.org for
> a complete listing of our services, staff and
> locations.
>
>
> Confidentiality Note:  This message is intended for use\...{{dropped:13}}
>
> __
> 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] political maps world maps in R, wrld_simpl

2009-04-06 Thread Patrizio Frederic
dear all,
I'm a newbie in map tools. I was asking to perform an apparently very
simple task.
I have a list of countries (about fifty all over in the world) scored
with a real positive value. Eg

Country score
Italy   .56
UK  .58
Korea   .41
Mexico  .63
...

I wish to plot geographical maps where every country is filled with a
color depending on his score.

Using maptools library and wrld_simpl I can easily perform the world
map. Here my code:

#
###plot the world map

patramp = colorRamp(c(rgb(0,0,.7),"red"))# my personal ramp color palette

library(maptools)
load(url("http://spatial.nhh.no/R/etc/TM_WORLD_BORDERS_SIMPL-0.2.RData";))

country2= c("Albania", "Argentina", "Australia", "Austria", "Belgium",
"Belarus", "Brazil", "Bulgaria", "Canada", "China", "Colombia",
"Croatia", "Denmark", "Egypt", "United.Arab.Emirates", "Estonia",
"Philippines", "Finland", "France", "Germany", "Japan", "Jordan",
"Greece", "India", "Indonesia", "Iran", "Ireland", "Israel",
"Italy", "Latvia", "Libya", "Lithuania", "Malaysia", "Mexico",
"Morocco", "Norway", "Netherlands", "Pakistan", "Poland", "Portugal",
"United.Kingdom", "Czech.Republic", "Romania", "Russia", "Serbia",
"Singapore", "Syrian.Arab.Republic", "Slovak.Republic", "Slovenia",
"Spain", "United.States", "South.Africa", "Korea.South", "Sweden",
"Switzerland", "Thailand", "Taiwan.Province.of.China", "Tunisia",
"Turkey", "Ukraine", "Hungary", "Venezuela", "Vietnam") # The my
countries' names

country.all = wrld_simpl$NAME # countries' names as coded in wrld_simpl
n.all   = length(country.all)
col.map = numeric(n.all)

# grep country2 in country.all
for (i in 1:length(country2)){
  col.map[grep(country2[i],country.all)] = runif(1) # assign a random
score to country i
}
# it certainly exists a cleaner coding for this but it's not my first
thought now

c2  = col.map!=0# set countries with score != 0 is not listed
in country2
col.map2=col.map
col.map2[!c2]=rgb(.95,.96,.97)  # be grey the countries not listed in
country2 object
col.map2[c2]= rgb(patramp(col.map[c2])/255) # ramp the remaining
plot(wrld_simpl,col=col.map2,axes=T) # nice plot

###   end of program
#

Now what if I wish to zoom on some sub world regions such as Eastern
Europe or America or any other user defined contiguous sub group?

Russia is a large country I could be interested in cutting it in two
peaces eastern and western. Can I do it using wrld_simpl or should I
search for other data? If so where can I find free maps data?

Some documentation exists but it seems sparse to me.
Any hints will be appreciated.
Thanks in advance,

Patrizio

__
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] political maps world maps in R, wrld_simpl

2009-04-07 Thread Patrizio Frederic
Roger,
many thanks for your very useful suggestions.
I've just produced some plot using xlim/ylim. I noticed I can only
zoom in a square box. Is there any way to produce rectangular ones?

Best

Patrizio

2009/4/7 Roger Bivand :
> Patrizio Frederic  gmail.com> writes:
>
>>
>> dear all,
>> I'm a newbie in map tools. I was asking to perform an apparently very
>> simple task.
>> I have a list of countries (about fifty all over in the world) scored
>> with a real positive value. Eg
>> ...
>>
>> I wish to plot geographical maps where every country is filled with a
>> color depending on his score.
>> Using maptools library and wrld_simpl I can easily perform the world
>> map.
>
> Note that the file is included in recent releases of maptools, so just say
>
> data(wrld_simpl)
>
>> library(maptools)
>> load(url("http://spatial.nhh.no/R/etc/TM_WORLD_BORDERS_SIMPL-0.2.RData";))
>
> For matching the countries you have with the countries in the object (note
> that names may differ), I suggest simply adding ISO3 codes to your data
> frame, and using them as row names as well. Look at the spCbind methods in
> maptools.
>
>> Now what if I wish to zoom on some sub world regions such as Eastern
>> Europe or America or any other user defined contiguous sub group?
>
> Either select using [ or plot using xlim/ylim arguments
>
>>
>> Russia is a large country I could be interested in cutting it in two
>> peaces eastern and western. Can I do it using wrld_simpl or should I
>> search for other data? If so where can I find free maps data?
>
> You can plot part of Russia with the object as it is by using xlim/ylim
> arguments. You can manipulate low-level objects if you need to, but
> it is not easy unless you know what you are doing.
>
>>
>> Some documentation exists but it seems sparse to me.
>
> The R-sig-geo list is more focussed on this kind of question, as is pointed
> out on the Spatial Task View on CRAN.
>
> Hope this helps,
>
> Roger Bivand
>
>> Any hints will be appreciated.
>> Thanks in advance,
>>
>> Patrizio
>
> __
> 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] Reversing axis label order

2009-05-04 Thread Patrizio Frederic
Dear Steve,
I'm not sure I properly understood the question. Try

x<- aggregate_1986 # just a shortcut
n<- dim(x)[[1]]
matplot(x[n:1,1], x[n:1,2:3], type="l", col=2:3)

is that what you needed?

Patrizio

2009/5/4 Steve Murray :
>
> Dear R Users,
>
> I am executing the following command to produce a line graph:
>
> matplot(aggregate_1986[,1], aggregate_1986[,2:3], type="l", col=2:3)
>
> On the x-axis I have values of Latitude (in column 1) ranging from -60 to +80 
> (left to right on the x-axis). However, I wish to have these values shown in 
> reverse on the x-axis, going from +80 to -60 (ie. North to South in terms of 
> Latitude). I have tried doing this by altering the command as follows:
>
> matplot(-aggregate_1986[,1], aggregate_1986[,2:3], type="l", col=2:3)
>
> ...but this produces the inverse sign of the latitude values along the axis - 
> ie. it goes from -80 to +60.
>
> How do I reverse the display of the axis labels correctly and of course, 
> maintain the associated data values correctly?
>
> Many thanks,
>
> Steve
>
> __
> 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.
>



-- 
+-
| Patrizio Frederic, PhD
| Assistant Professor,
| Department of Economics,
| University of Modena and Reggio Emilia,
| Via Berengario 51,
| 41100 Modena, Italy
|
| tel:  +39 059 205 6727
| fax:  +39 059 205 6947
| mail: patrizio.frede...@unimore.it
+-

__
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 create a numeric data.frame

2011-06-13 Thread Patrizio Frederic
On Mon, Jun 13, 2011 at 4:45 PM, Barry Rowlingson
 wrote:
> On Mon, Jun 13, 2011 at 11:06 AM, Aparna  wrote:
>> Hi All
>>
>> I am new to R and  I am not sure of how this should be done. I have a matrix 
>> of
>> 985x100 values and the class is data.frame.
>
>  You don't have a 'matrix' in the R sense of the word. You seem to
> have a table of numbers which are stored in an object of class
> 'data.frame'.
>

but you could have one:

a <- data.frame(matrix(rnorm(100),10) # get some data
class(a) # check for its class
as.numeric(a) # whoops, won't work

class(as.matrix(a)) # change class, and
as.numeric(as.matrix(a)) # bingo, it works

PF


-- 
+---
| Patrizio Frederic,
| http://www.economia.unimore.it/frederic_patrizio/
+---

__
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 create a numeric data.frame

2011-06-13 Thread Patrizio Frederic
>
> Which results in vector of numbers
>
> str(as.numeric(as.matrix(a)))
>  num [1:100] 0.82 -1.339 1.397 0.673 -0.461 ...
>
> data frame is convenient list structure which can contain vectors of
> various nature (numeric, character, factor, logical, ...)
> and looks quite similar to Excel table.
>
> matrix is a vector with (2) dimensions but as it is a vector it can not
> consist from objects of different nature (class). Therefore you can have
> numeric or character matrix but not numeric and character columns in your
> matrix.
>
> and vector is vector (numeric, character, logical,  ...) but again you can
> not mix items of different class in one vector.
>

of course it is. I forgot to say that the way I proposed works only if
the data-frame contains numeric objects only.

R is a great tool because you can get to the very same results in many
different ways.
Depending on the problem you're dealing with, you have to choose the
most efficient one.
Often, in my research work, the most efficient is the one that use as
less as possible lines of code:

Suppose a is a data.frame which contains numeric objects only

a <- data.frame(matrix(rnorm(100),10)) # some data

## 1 not very nice
b <- 0
for (j in 1:length(a)) b<-c(b,as.numeric(a[i]))
b<-b[-1]

## 2 long time ago I was a fortran guy
b<-numeric(length(a))
for (j in 1:dim(a)[2]){
  for (i in 1:dim(a)[1]){
 b[10*(j-1)+i] <- as.numeric(a[i,j])
  }
}

## 3 better: sapply function
as.numeric(sapply(a,function(x)as.numeric(x)))

## 4 shorter
as.numeric(as.matrix(a))

## which type of data a has
a <- data.frame(a,fact=sample(c('F1','F2'),dim(a)[1],replace=T))
class_a <- sapply(a,function(x)class(x))
class_a
a_numeric <- a[,class_a=='numeric']
as.numeric(as.matrix(a_numeric))

Regards,

PF

-- 
+---
| Patrizio Frederic,
| http://www.economia.unimore.it/frederic_patrizio/
+---

__
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 create a numeric data.frame

2011-06-13 Thread Patrizio Frederic
On Mon, Jun 13, 2011 at 6:47 PM, Aparna  wrote:
> Hi Joshua
>
> While looking at the data, all the values seem to be in numeric. As i 
> mentioned,
> the dataset is already in data.frame.
>
> As suggested, I used str(mydata) and got the following result:
>
>
> str(leu_cluster1)
> 'data.frame':   984 obs. of  100 variables:
>  $ V2  : Factor w/ 986 levels "-0.00257361",..: 543 116 252 54 520 ...

your data columns are not numeric but factors indeed.
you may try this one

a <- as.character(rnorm(100))   # some numeric data
adf <- data.frame(matrix(a,10)) # which are misinterpreted as factors
adf
adf[,1]
class(adf[,1]) # check for the class of the first column
sapply(adf,function(x)class(x)) # check classes for all columns

b <- sapply(adf,function(x)as.numeric(as.character(x))) #
as.character: use levels literally, as.numeric: transforms in numbers
b # look at b

class(b) # which is now a numeric matrix

best regards

PF

-- 
+---
| Patrizio Frederic,
| http://www.economia.unimore.it/frederic_patrizio/
+---

__
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] question about R

2010-05-11 Thread Patrizio Frederic
# density function
desp <- function(x,lambda) {
  return(lambda * exp(-lambda*x)*(x>=0))
}
# cum. probability function
pesp <- function(q,lambda) {
  return((1-exp(-lambda*q)))
}
# quantile function
qesp <- function(p,lambda) {
  return(-lambda^{-1}*log(1-p))
}
# random sample
resp <- function(n,lambda) {
  return(qesp(runif(n),lambda))
}
resp(100,lambda=2) # generates 100 samples from the exponential
distribution with par 2


On Tue, May 11, 2010 at 6:36 PM,   wrote:
> Hi,
> At each iteration in my program,I need to generate  tree vectors,X1,X2,X3,
> from exponential distribution with parameters a1,a2,a3. Can you help me
> please how can I do it such that it take a little time?
>
> thank you
> khazaei
>
> __
> 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.
>



-- 
+-----
| Patrizio Frederic, PhD
| Assistant Professor,
| Department of Economics,
| University of Modena and Reggio Emilia,
| Via Berengario 51,
| 41100 Modena, Italy
|
| tel:  +39 059 205 6727
| fax:  +39 059 205 6947
| mail: patrizio.frede...@unimore.it
+-

__
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] access objects in my environment

2010-05-13 Thread Patrizio Frederic
Hi,
have a look at function 'get'  (and function 'assign' too)
In your example:

x <- 1000
assign(paste('Pos',x,sep=''),rnorm(100))
Pos1000
posa <- get(paste(c("Pos",x),collapse=""))
posa

hope it's help

PF

On Thu, May 13, 2010 at 10:37 AM, arnaud Gaboury
 wrote:
> Dear group,
>
> Here are my objects in my environment:
>
>> ls()
> [1] "Pos100415" "Pos100416" "posA"      "pose15"    "pose16"    "pose16t"
> "position"  "trade"     "x"
>
> I need to pass the object "Pos100415" to a function. This element is a
> data.frame, obtained through a function: Pos(x)<-myfun(x) with x<-100415 in
> this example.
>
> If I do this :
>
>>posA<-paste(c("Pos",100415),collapse="") #I can't use directly Pos100415.I
> need to access it via paste(c("Pos",x),collapse="") in general.
>
> Here is what I got :
>
>>> posA
> [1] "Pos100415"
>
> It is certainly not what I want, as I need to have posA as the same
> data.frame than Pos100415.
>
> Any help?
>
> TY
>
> ______
> 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.
>



-- 
+-
| Patrizio Frederic, PhD
| Assistant Professor,
| Department of Economics,
| University of Modena and Reggio Emilia,
| Via Berengario 51,
| 41100 Modena, Italy
|
| tel:  +39 059 205 6727
| fax:  +39 059 205 6947
| mail: patrizio.frede...@unimore.it
+-

__
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] best polynomial approximation

2010-05-17 Thread Patrizio Frederic
Dear R-users,
I learned today that there exists an interesting topic in numerical
analysis names "best polynomial approximation" (BSA). Given a function
f  the BSA of degree k, say pk, is the polynomial such that

pk=arginf sup(|f-pk|)

Although given some regularity condition of f, pk is unique, pk IS NOT
calculated with least square. A quick google tour show a rich field of
research and many algorithms proposed for computing such a task.

I was wondered if some of you knows about some R implementations
(packages) for computing BSA.

Many thanks in advance,

Patrizio

as usual I apologize for my fragmented English

-- 
+-----
| Patrizio Frederic, PhD
| Assistant Professor,
| Department of Economics,
| University of Modena and Reggio Emilia,
| Via Berengario 51,
| 41100 Modena, Italy
|
| tel:  +39 059 205 6727
| fax:  +39 059 205 6947
| mail: patrizio.frede...@unimore.it
+-

__
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] best polynomial approximation

2010-05-18 Thread Patrizio Frederic
Dear colleagues,
thank you so much for your help.
Hans, I think the Remez algorithm is what I need. I will brush up on
fortran language.
Ravi, thanks anyway, I appreciated.

All the best,

Patrizio



On Tue, May 18, 2010 at 12:10 PM, Hans W Borchers
 wrote:
>
> I guess you may be looking for the Remez algorithm. AFAIK there is no
> implementation in one of the R packages. You can find FORTRAN code in the
> Collected Algorithms of the ACM (no. 604) which probably could be called
> from R.
>
> There appears to exist a discrete, equi-distant(?) version as function
> 'remez' in the signal package, if that is of any help to you. I have never
> used it.
>
> Regards,  Hans Werner
>
> P.S.: The Chebyshev polynomials do not compute the "best polynomial
> approximation", but they provide a nice way to estimate the maximal distance
> to this best approximating polynomial.
>
>
>
> Patrizio Frederic wrote:
>>
>> Dear R-users,
>> I learned today that there exists an interesting topic in numerical
>> analysis names "best polynomial approximation" (BSA). Given a function
>> f  the BSA of degree k, say pk, is the polynomial such that
>>
>> pk=arginf sup(|f-pk|)
>>
>> Although given some regularity condition of f, pk is unique, pk IS NOT
>> calculated with least square. A quick google tour show a rich field of
>> research and many algorithms proposed for computing such a task.
>>
>> I was wondered if some of you knows about some R implementations
>> (packages) for computing BSA.
>>
>> Many thanks in advance,
>>
>> Patrizio
>>
>> as usual I apologize for my fragmented English
>>
>> --
>> +-
>> | Patrizio Frederic, PhD
>> | Assistant Professor,
>> | Department of Economics,
>> | University of Modena and Reggio Emilia,
>> | Via Berengario 51,
>> | 41100 Modena, Italy
>> |
>> | tel:  +39 059 205 6727
>> | fax:  +39 059 205 6947
>> | mail: patrizio.frede...@unimore.it
>> +-
>>
>> __
>> 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.
>>
>>
> --
> View this message in context: 
> http://r.789695.n4.nabble.com/best-polynomial-approximation-tp2220439p2221042.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.
>



-- 
+-
| Patrizio Frederic, PhD
| Assistant Professor,
| Department of Economics,
| University of Modena and Reggio Emilia,
| Via Berengario 51,
| 41100 Modena, Italy
|
| tel:  +39 059 205 6727
| fax:  +39 059 205 6947
| mail: patrizio.frede...@unimore.it
+-

__
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] Finding non-normal distributions per row of data frame?

2011-02-04 Thread Patrizio Frederic
Hi Danny,
it sounds to much easer than that.

Try

y <- data.frame(matrix(rnorm(100),10))
nr <- ncol(y)

test <- lapply(y, shapiro.test)
sapply(test,function(x)c(x$statistic, x$p.value))

it should perform the required task.
Cheers,

P

On Fri, Feb 4, 2011 at 4:52 AM, DB1984  wrote:
>
> This is my first attempt at this, so hopefully a few kind pointers can get me
> going in the right direction...
>
> I have a large data frame of 20+ columns and 20,000 rows. I'd like to
> evaluate the distribution of values in each row, to determine whether they
> meet the criteria of a normal distribution. I'd loop this over all the rows
> in the data frame, and output the summary results to a new data frame.
>
> I have a loop that should run a Shapiro-Wilk test over each row,
>
> y= data frame
>
> for (j in 1:nr) {
> y.temp<-list(y[j,])
> testsw <- lapply(y.temp, shapiro.test)
> testtable <- t(sapply(testsw, function(x) c(x$statistic, x$p.value)))
>  colnames(testtable) <- c("W", "p.value")
> }
>
>
> but it is currently throwing out an error:
>  "Error in `rownames<-`(`*tmp*`, value = "1") :
>  attempt to set rownames on object with no dimensions"
>
> ...which I guess is unrelated to the evaluation of normality, and more
> likely a faulty loop?
>
> Any suggestions either for this test, or a better way to evaluate the normal
> distribution (e.g. qq-plot residuals for each row) would be greatly
> received. Thanks.
> --
> View this message in context: 
> http://r.789695.n4.nabble.com/Finding-non-normal-distributions-per-row-of-data-frame-tp3259439p3259439.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] normal distribution and floating point traps (?): unexpected behavior

2011-03-29 Thread Patrizio Frederic
dear all,
here's a couple of questions that puzzled me in these last hours:

# issue 1 qnorm(1-10e-100)!=qnorm(10e-100)

qnorm(1-1e-10) == -qnorm(1e-10)

# turns on to be FALSE. Ok I'm not a computer scientist but,
# but I had a look at the R inferno so I write:

all.equal(qnorm(1-1e-10) , -qnorm(1e-10))

# which turns TRUE, as one would expect, but

all.equal(qnorm(1-1e-100) , -qnorm(1e-100))

# turns FALSE: Indeed

# qnorm(1e-100) is -21.27345, and
# qnorm(1-1e-100) is Inf

# as a consequence there was an infinitive (I would expect a 21.27)
running in my program which was very annoying and hard to detect.

# issue 2 dnorm(...,log=logical_flag)
# simple normal likelihood

set.seed(1)
x <- rnorm(100) # sample 100 normal data
mu <- seq(-4,4,length=51) # get a grid for mu
sigma <- seq(.01,4,length=51) # get a grid for sigma
lik <- function(theta,x) sum( dnorm(x,theta[1],theta[2],log=T)) # log likelihood
Lik <- function(theta,x) prod(dnorm(x,theta[1],theta[2],log=F)) # Likelihood
ms  <- as.matrix(expand.grid(mu,sigma))

l   <- matrix(apply(ms,1,lik,x=x),51)
L   <- matrix(apply(ms,1,Lik,x=x),51)

contour(mu,sigma,L) # plots exactly what I expected
contour(mu,sigma,log(L)) # plots exactly what I expected
contour(mu,sigma,l) # I didn't expect that!!!'
contour(mu,sigma,log(exp(l))) # plots exactly what I expected

> version
   _
platform   x86_64-pc-linux-gnu
arch   x86_64
os linux-gnu
system x86_64, linux-gnu
status
major  2
minor  12.2
year   2011
month  02
day25
svn rev54585
language   R
version.string R version 2.12.2 (2011-02-25)

under ubuntu distribution

My questions are: is that a bug? Does this behavior depend on my pc
architecture?
What kind of precautions does R's guru suggests to be taken?

Thank you in advance,

PF

ps as usual please I apologize for my fragemnted English

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


Re: [R] R PLS package data format

2011-06-27 Thread Patrizio Frederic
hi,
because of my ignorance I can't figure out why would you need such a structure.
However the following trick works:

> a <- matrix(1:6,3)
> b <- data.frame(1:3)
> b$a <- a
> b
  X1.3 a.1 a.2
11   1   4
22   2   5
33   3   6
> names(b)
[1] "X1.3" "a"
> class(b$a)
[1] "matrix"
> is.matrix(b$a)
[1] TRUE

hope that's help

PF

2011/6/27 新鼎-智慧製造事業處-鄭紹文 :
> Dear sir,
>
> If I have a vector as:
>>y <- c(1:4)
> and a matrix as:
>>x <- matrix (5:12, nrow=4, ncol=2)
>
> Then I create a data frame as:
>>t <- data.frame(y, x)
>
> R will generate a data frame of 3 columns
>>names(t)
>> "y"  "X1" "X2"
>>is.vector(t$X1) returns TRUE
>> is.vector(t$X2) returns TRUE
> R splits the columns of a matrix into column vectors automatically.
>
> How should I create a data frame of t$y & t$x such that
>>is.vector(t$y) returns TRUE
>>is.matrix(t$x) returns TRUE
> That is, how to make x column of the data frame a matrix? Just like the test 
> data ¨gasoline〃 in your PLS package:
>>library(pls)
>>data(yarn)
>>is.vector(density) returns TRUE
>>is.matrix(NIR) returns TURE (NIR is a matrix)
>
> Best regards,
> Shaowen Cheng
>
> 
> The information in this e-mail may be confidential; it is intended for use 
> solely by the individual or entity named as the recipient hereof. Disclosure, 
> copying, distribution, or use of the contents of this e-mail by persons other 
> than the intended recipient may violate applicable laws and if you have 
> received this e-mail in error, please delete the original message and notify 
> us by collect call immediately. Thank you.
>
>        [[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.
>
>



-- 
+---
| Patrizio Frederic,
| http://www.economia.unimore.it/frederic_patrizio/
+---

__
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] Improve a browse through list items - Transform a loop to apply-able function.

2011-12-14 Thread Patrizio Frederic
Hi robin,
I'm not sure is what you need, but that's an esthetically nice
solution (one single line without any loop :) )

matrix(apply(log(cbind(as.numeric(a),as.numeric(b),as.numeric(c),as.numeric(d))),1,sd),3)

hope it could help,

PF


On Mon, Dec 12, 2011 at 5:15 PM, Robin Cura  wrote:
> Hello,
>
> I'm currently trying to convert a slow and ugly script I made, so that it's
> faster and can be computed on a computer grid with the multicore package.
> My problem is that I don't see how to turn some loops into an "apply-able"
> function.
>
> Here's an example of my loops :
> I got a list of dataframes (or matrices like here), and I need to browse
> each cell of those many dataframes to compute a mean (or standard deviation
> like here).
>
> Here's a example script :
>
> a <- b <- c <- d <- result <- matrix(nrow=3, ncol=3)
> a[] <- sample.int(n=100,size=9,replace=TRUE)
> b[] <- sample.int(n=100,size=9,replace=TRUE)
> c[] <- sample.int(n=100,size=9,replace=TRUE)
> d[] <- sample.int(n=100,size=9,replace=TRUE)
> result[] <- NA
> mylist <- list(a,b,c,d)
>
> for (row in 1:3)
> {
>  for (col in 1:3)
>  {
>    tmpList <- log(mylist[[1]][row, col])
>    for (listitem in 2:4)
>    {
>      tmpList <- c(tmpList, log(mylist[[listitem]][row, col]))
>    }
>    result[row, col] <- sd(tmpList)
>  }
> }
>
> Considering I have to look at the same cell in each dataframe, I don't
> understand how I could turn this into a function, considering I need the
> row and column number to iterate.
>
> I succeeded improving my script duration a lot, but such loops are really
> long to run, considering that my lists contains like 100 dataframes, who
> all contains thousands of values.
>
> Any help would be really appreciated
>
> Thanks in advance,
>
> Robin
>
>        [[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.



-- 
+---
| Patrizio Frederic,
| http://www.economia.unimore.it/frederic_patrizio/
+---

__
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] While loop working with TRUE/FALSE?

2012-02-02 Thread Patrizio Frederic
Hey Chris,
I would take advantage from the apply function:

apply(cbind(q1,q2),1,function(x)any((x[1]==w1)&(x[2]==w2)))

Regards

PF

On Thu, Feb 2, 2012 at 12:55 PM, Chris82  wrote:
> Thanks to Berend and the others,
>
> I've found a solution which works fine for my problem.
>
> I have not only 2 vectors, but also 4.
> Question is, if q1 and q2 is equal to w1 and w2.
> The computational time is very short, also for large data.
>
> q1 <- c(9,5,1,5)
> q2 <- c(9,2,1,5)
>
> w1 <- c(9,4,4,4,5)
> w1 <- c(9,4,4,4,5)
>
> v <- vector()
> for (i in 1:(length(q1))){
> v[i] <- any((q1[i] == w1) &  (q2[i] == w2))
> }
>
>
> best regards
>
> --
> View this message in context: 
> http://r.789695.n4.nabble.com/While-loop-working-with-TRUE-FALSE-tp4348340p4351214.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.



-- 
+---
| Patrizio Frederic,
| http://www.economia.unimore.it/frederic_patrizio/
+---

__
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] Troubles plotting lrm output in Design Library

2008-05-29 Thread Patrizio Frederic
Dear R-helpers,
I'm having a problem in using plot.design in Design Library. Tho
following example code produce the error:

> n <- 1000# define sample size
> set.seed(17) # so can reproduce the results
> age<- rnorm(n, 50, 10)
> blood.pressure <- rnorm(n, 120, 15)
> cholesterol<- rnorm(n, 200, 25)
> sex<- factor(sample(c('female','male'), n,TRUE))
> label(age)<- 'Age'  # label is in Hmisc
> label(cholesterol)<- 'Total Cholesterol'
> label(blood.pressure) <- 'Systolic Blood Pressure'
> label(sex)<- 'Sex'
> units(cholesterol)<- 'mg/dl'   # uses units.default in Hmisc
> units(blood.pressure) <- 'mmHg'
>
> # Specify population model for log odds that Y=1
> L <- .4*(sex=='male') + .045*(age-50) +
+   (log(cholesterol - 10)-5.2)*(-2*(sex=='female') + 2*(sex=='male'))
> # Simulate binary y to have Prob(y=1) = 1/[1+exp(-L)]
> y <- ifelse(runif(n) < plogis(L), 1, 0)
>
> ddist <- datadist(age, blood.pressure, cholesterol, sex)
> options(datadist='ddist')
>
> fit <- lrm(y ~ blood.pressure + sex * (age + rcs(cholesterol,4)),
+x=TRUE, y=TRUE)
>
> par(mfrow=c(2,2))
> plot(fit)# Plot effects of all 4 predictors
Error in value.chk(at, ix, xseq, if (plot.type == "curves") 100 else 40,  :
 variable blood.pressure does not have limits defined by datadist

I guess somthing goes wrong with function datadist

> version
  _
platform   i386-pc-mingw32
arch   i386
os     mingw32
system i386, mingw32
status
major      2
minor      7.0
year   2008
month  04
day22
svn rev45424
language   R
version.string R version 2.7.0 (2008-04-22)

thank you in advance,

Patrizio Frederic

-
| Patrizio Frederic
| Research associate in Statistics,
| Department of Economics,
| University of Modena and Reggio Emilia,
| Via Berengario 51,
| 41100 Modena, Italy
|
| tel:  +39 059 205 6727
| fax:  +39 059 205 6947
| mail: [EMAIL PROTECTED]
+-

__
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] Troubles plotting lrm output in Design Library

2008-05-29 Thread Patrizio Frederic
dear Harrell,
thank you for quick reply and suggestions. I still have the problem:

library(Design)
x   = rnorm(100)
y   = runif(100)<(exp(x)/(1+exp(x)))
y   = 0*y+1*y
d   = datadist(x,y)
options(datadist="d")
fit = lrm(y~x)
# works fine, but
plot(fit) #produce the error message

Error in value.chk(at, ix, xseq, if (plot.type == "curves") 100 else 40,  :
  variable x does not have limits defined by datadist

The installed Design (ver. 2.1.1) package was downloaded and installed
via utils:::menuInstallPkgs() function
(Package Hmisc version 3.4-3) on a windows xp machine.
Thank in advance.

Patrizio Frederic

version

platform   i386-pc-mingw32
arch   i386
os mingw32
system i386, mingw32
status
major  2
minor  7.0
year   2008
month  04
day22
svn rev45424
language   R
version.string R version 2.7.0 (2008-04-22)

2008/5/29 Frank E Harrell Jr <[EMAIL PROTECTED]>:
> Patrizio Frederic wrote:
>>
>> Dear R-helpers,
>> I'm having a problem in using plot.design in Design Library. Tho
>> following example code produce the error:
>>
>>> n <- 1000# define sample size
>>> set.seed(17) # so can reproduce the results
>>> age<- rnorm(n, 50, 10)
>>> blood.pressure <- rnorm(n, 120, 15)
>>> cholesterol<- rnorm(n, 200, 25)
>>> sex<- factor(sample(c('female','male'), n,TRUE))
>>> label(age)<- 'Age'  # label is in Hmisc
>>> label(cholesterol)<- 'Total Cholesterol'
>>> label(blood.pressure) <- 'Systolic Blood Pressure'
>>> label(sex)<- 'Sex'
>>> units(cholesterol)<- 'mg/dl'   # uses units.default in Hmisc
>>> units(blood.pressure) <- 'mmHg'
>>>
>>> # Specify population model for log odds that Y=1
>>> L <- .4*(sex=='male') + .045*(age-50) +
>>
>> +   (log(cholesterol - 10)-5.2)*(-2*(sex=='female') + 2*(sex=='male'))
>>>
>>> # Simulate binary y to have Prob(y=1) = 1/[1+exp(-L)]
>>> y <- ifelse(runif(n) < plogis(L), 1, 0)
>>>
>>> ddist <- datadist(age, blood.pressure, cholesterol, sex)
>>> options(datadist='ddist')
>>>
>>> fit <- lrm(y ~ blood.pressure + sex * (age + rcs(cholesterol,4)),
>>
>> +x=TRUE, y=TRUE)
>>>
>>> par(mfrow=c(2,2))
>>> plot(fit)# Plot effects of all 4 predictors
>>
>> Error in value.chk(at, ix, xseq, if (plot.type == "curves") 100 else 40,
>>  :
>>  variable blood.pressure does not have limits defined by datadist
>>
>> I guess somthing goes wrong with function datadist
>>
>>> version
>>
>>  _
>> platform   i386-pc-mingw32
>> arch   i386
>> os mingw32
>> system i386, mingw32
>> status
>> major  2
>> minor  7.0
>> year   2008
>> month  04
>> day22
>> svn rev45424
>> language   R
>> version.string R version 2.7.0 (2008-04-22)
>>
>> thank you in advance,
>>
>> Patrizio Frederic
>
> This is the example run by example(plot.Design) which works for me, using
> the version of Design from CRAN.  You don't need to include code that is
> already in an example in a help file, and if you do please make the code
> copy and paste-able instead of putting something at the start of each line.
>
> Frank
>
>
> Frank
>
>>
>> -
>> | Patrizio Frederic
>> | Research associate in Statistics,
>> | Department of Economics,
>> | University of Modena and Reggio Emilia,
>> | Via Berengario 51,
>> | 41100 Modena, Italy
>> |
>> | tel:  +39 059 205 6727
>> | fax:  +39 059 205 6947
>> | mail: [EMAIL PROTECTED]
>> +-
>>
>> __
>> 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.
>>
>
>
> --
> Frank E Harrell Jr   Professor and Chair   School of Medicine
> Department of Biostatistics   Vanderbilt University
>

__
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] Troubles plotting lrm output in Design Library

2008-05-30 Thread Patrizio Frederic
2008/5/29 Frank E Harrell Jr <[EMAIL PROTECTED]>:
> Patrizio Frederic wrote:
>>
>> dear Harrell,
>> thank you for quick reply and suggestions. I still have the problem:
>>
>> library(Design)
>> x   = rnorm(100)
>> y   = runif(100)<(exp(x)/(1+exp(x)))
>> y   = 0*y+1*y
>> d   = datadist(x,y)
>> options(datadist="d")
>> fit = lrm(y~x)
>> # works fine, but
>> plot(fit) #produce the error message
>
> I cannot reproduce the error on Linux R 2.7.0 using the same version of the
> packages you are using.
>
> Frank
>

just for the records, I solved the problem: in Windows R 2.7.0 one has to set:

d = datadist(x,y)
.Options$datadist="d"

done that functions summary.Design, and plot.Design work fine.

PF

__
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] Troubles plotting lrm output in Design Library

2008-05-30 Thread Patrizio Frederic
> Thanks for the solution.  But I wish someone would tell us why the solution
> makes sense.  If you have changed your R environment in any way from the
> defaults please let us know.
>
> Frank
>
> --
> Frank E Harrell Jr   Professor and Chair   School of Medicine
> Department of Biostatistics   Vanderbilt University
>

dear Harrell and R users,
It took me several hours to figure out. I apologize it was my mistake.
Few month ago a modified copy of .Options was created in Rprofile (I
don't remember why). The copy masked the alias .Options created by
options().
Now everything is fixed. Now I know how .Options and options are
related one each other. Sorry for wasting your time.
Thank you for your help.

Patrizio Frederic

__
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] migrating data from s-plus to R

2008-09-19 Thread Patrizio Frederic
Dear all,
is there any way to transform a .Data directory created in S-plus 6.1
for windows in a .RData file?
Thanks in advance,

Patrizio Frederic

__
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] multiple secondary axes

2009-01-14 Thread Patrizio Frederic
hey Kirsten,
did type

?axis

it may help. eg:

plot(rnorm(100),axes=F)
axis(1)
axis(3,at=seq(0,50,length=3),labels=c("A","B","C"))
axis(2)
axis(4,at=seq(-2,2,length=10),labels=1:10)

Best,

Patrizio

2009/1/14 Kirsten Thonicke :
> Dear R experts,
>   I want to plot a line chart with another secondary axis placed right to
> the standard secondary axis which one can access with the axis command, so
> that the data lines are seen in the same plot. Is there any way to do this
> in R?
> Many thanks,
> Kirsten.
>
> __
> 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] extracting parts of words or extraxting letter to use in ifelse-func.

2009-02-13 Thread Patrizio Frederic
hey Joe,
I had a quick look at your code. In this line:

b2$PRRSvac <- ifelse(b2$status=="MS-X, PRRS-pos"|b2$status=="R�d
SPF+Myc+Ap2+DK+Vac",1,0)

some special characters are used. You must use only plain ascii characters
Hope this help.

Cheers,

Patrizio

2009/2/13 joe1985 :
>
> Hello
>
> I want to make some variables with the ifelse-function, but i don't know how
> to do it.
>
> I want to make these five variables;
>
> b2$PRRSvac <- ifelse(b2$status=='A' | b2$status=='Aa',1,0)
> b2$PRRSdk <- ifelse(b2$status=='B' | b2$status=='Bb',1,0)
> b2$sanVac <- ifelse(b2$status=='C' | b2$status=='sanAa',1,0)
> b2$sanDk <- ifelse(b2$status=='D'  | b2$status=='sanBb',1,0)
> b2$totalvac <- ifelse(b2$status=='San',1,0)
>
> And
>
> b2$UGT <- ifelse(b2$status=='UGT',1,0)
> b2$KOM <- ifelse(b2$status=='KOM',1,0)
>
> But, as one from this forum told me, it doesn't work because the words is in
> a wrong format or something like that.
>
> I have attached the text-file i've used, and the R-kode.
>
> Hope anyone can help me?
>
> text-file;
> http://www.nabble.com/file/p21993098/allesd%2528uden%2Bf%25C3%25B8r%2529060209.txt
> allesd%28uden+f%C3%B8r%29060209.txt
>
> R-kode (just run from line 1-22) ;
> http://www.nabble.com/file/p21993098/PRRS%2B%2528med%2Ballle%2BSD%2529.r
> PRRS+%28med+allle+SD%29.r
> --
> View this message in context: 
> http://www.nabble.com/extracting-parts-of-words-or-extraxting-letter-to-use-in-ifelse-func.-tp21993098p21993098.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] several "ifelse" problems...

2009-02-16 Thread Patrizio Frederic
do you mean:

 f=function(x) 
0*(abs(x-.5)<=.3)-1*(abs(x-.5)>=.4)+(10*x-2)*(x>.1&x<.2)+(-10*x+8)*(x>=.2&x<=.5)
 f(x)
 curve(f,0,1)

hope it helps.

Patrizio

2009/2/14 kathie :
>
> Dear R users,
>
> >From the code below, I try to compute "y" value. (In fact, y looks like a
> trapezoid)
>
> --
>
> x <- seq(0,1,.01)
> y <- ifelse(abs(x-.5)<=0.3,0,
>ifelse(abs(w-.5)>=0.4,-1,
>   ifelse((0.1
> --
>
> So, results are...
>
> --
>> x
>  [1] 0.00 0.01 0.02 0.03 0.04 0.05 0.06 0.07 0.08 0.09 0.10 0.11 0.12 0.13
> 0.14
>  [16] 0.15 0.16 0.17 0.18 0.19 0.20 0.21 0.22 0.23 0.24 0.25 0.26 0.27 0.28
> 0.29
>  [31] 0.30 0.31 0.32 0.33 0.34 0.35 0.36 0.37 0.38 0.39 0.40 0.41 0.42 0.43
> 0.44
>  [46] 0.45 0.46 0.47 0.48 0.49 0.50 0.51 0.52 0.53 0.54 0.55 0.56 0.57 0.58
> 0.59
>  [61] 0.60 0.61 0.62 0.63 0.64 0.65 0.66 0.67 0.68 0.69 0.70 0.71 0.72 0.73
> 0.74
>  [76] 0.75 0.76 0.77 0.78 0.79 0.80 0.81 0.82 0.83 0.84 0.85 0.86 0.87 0.88
> 0.89
>  [91] 0.90 0.91 0.92 0.93 0.94 0.95 0.96 0.97 0.98 0.99 1.00
>
>> y
>  [1] -1 -1 -1 -1 -1 -1 -1 -1 -1 -1 -1  8  8  8  8  8  8  8  8  8  0  0  0
> 0  0
>  [26]  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0
> 0  0
>  [51]  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0
> 0  0
>  [76]  0  0  0  0  0  8  8  8  8  8  8  8  8  8  8 -1 -1 -1 -1 -1 -1 -1 -1
> -1 -1
> [101] -1
>>
>
> --
>
> However, even though the results show that y=8 for x=0.11, when x=0.11,
> actual y value is -0.9.  And, y=-0.8 for x=0.88.  I cannot understand the
> above results.
>
> Any comments will be greatly appreciated.
>
> Kathryn Lord
> --
> View this message in context: 
> http://www.nabble.com/several-%22ifelse%22-problems...-tp22009321p22009321.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] Unadulterated plot

2009-02-18 Thread Patrizio Frederic
James,
you're probably interested in image function rather than in
filled.contour. Type ?image to see the syntax.
Cheers,

Patrizio

2009/2/18 James Nicolson :
> Hi,
>
> Thanks for your help. I have looked at the beginners documentation and
> while there are options to configure various aspects of the plot none of
> them seem to have the desired effect. I have managed to ensure that the
> plot fills the space vertically with no margins, no axes etc (using
> mai=c(0,0,0,0)). However, horizontally there remains a margin to the
> right that pads the space between the filled.contour and its legend.
> I've tried options to par and filled.contour but I can't seem to remove
> the legend.
>
> Kind Regards,
> James
>
> Simon Pickett wrote:
>> Hi James,
>>
>> What you really need to do is to check out the many freely available
>> pdfs for R beginners. Here is a good place to start
>>
>> http://cran.r-project.org/other-docs.html
>>
>> If I am right interpreting what you want, I think you need to create a
>> blank plot with no axes, axis labels etc. Try
>>
>> plot(x,y,xlab="",ylab="",xaxt=NULL,yaxt=NULL,type="n")
>> #blank plot
>> points(x,y)
>>
>> type "?par" into R and see how you can set parameters like this up as
>> the default.
>>
>> Hope this helps?
>>
>> Simon.
>>
>>
>> - Original Message - From: "James Nicolson"
>> 
>> To: 
>> Sent: Sunday, February 15, 2009 10:29 PM
>> Subject: [R] Unadulterated plot
>>
>>
>>> To all,
>>>
>>> Apologies if this question has already been asked but I can't find
>>> anything. I can't seem to think of more specific search terms. I want
>>> to display/create a file of a pure plot with a specific height and
>>> width. I want to utilise every single pixel inside the axes. I do not
>>> want to display any margins, legends, axes, titles or spaces around
>>> the edges. Is this possible? Additionally, the plot I am working with
>>> is a filled.contour plot and I can not remove the legend? How can I
>>> do this?
>>>
>>> Kind Regards,
>>> James
>>>
>>> __
>>> 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-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] Unadulterated plot

2009-02-19 Thread Patrizio Frederic
James,
as I previously told you in my broken  English, probably the function
you're looking for is not filled.contour but image and contour
The following code makes exactly what you ask for

data(akima)
akima
akima.smooth <-
 with(akima, interp(x, y, z, xo=seq(0,25, length=500),
yo=seq(0,20, length=500)))
op <- par(ann=FALSE, mai=c(0,0,0,0))
image  (akima.smooth, main = "interp(, *) on finer grid")
contour(akima.smooth, add = TRUE,drawlabels=F)

cheers

Patrizio

2009/2/19 James Nicolson :
> good point! Provide your own set of x,y,z co-ords, mine are pretty big
> but you can use any.
>
> library(akima)
>
> fr3d = data.frame(x,y,z)
> xtrp <- interp(fr3d$x,fr3d$y,fr3d$z,linear=FALSE,extrap=TRUE,duplicate=
> "strip")
>
> op <- par(ann=FALSE, mai=c(0,0,0,0))
> filled.contour(xtrp$x, xtrp$y, xtrp$z, asp = 0.88402366864, col =
> rev(rainbow(28,start=0, end=8/12)), n = 40)
> par(op)
>
> I tried all these settings too (none of them made a difference)...
> usr=c(0,845,0,747), mfcol=c(1,1), mfrow=c(1,1),
> oma=c(0,0,0,0),omi=c(0,0,0,0), plt=c(1,1,1,1)
>
> Regards
> James
>
>
> Peter Alspach wrote:
>> Kia ora James
>>
>> I think it would be easier to provide you with help if you "provide
>> commented, minimal, self-contained, reproducible code" [see bottom of
>> this, or any, email to R-help].
>>
>> Hei kona ra ...
>>
>> Peter Alspach
>>
>>
>>> -Original Message-
>>> From: r-help-boun...@r-project.org
>>> [mailto:r-help-boun...@r-project.org] On Behalf Of James Nicolson
>>> Sent: Thursday, 19 February 2009 11:22 a.m.
>>> To: r-help@r-project.org
>>> Subject: Re: [R] Unadulterated plot
>>>
>>> Hi,
>>>
>>> Thanks for your help. I have looked at the beginners
>>> documentation and while there are options to configure
>>> various aspects of the plot none of them seem to have the
>>> desired effect. I have managed to ensure that the plot fills
>>> the space vertically with no margins, no axes etc (using
>>> mai=c(0,0,0,0)). However, horizontally there remains a margin
>>> to the right that pads the space between the filled.contour
>>> and its legend.
>>> I've tried options to par and filled.contour but I can't seem
>>> to remove the legend.
>>>
>>> Kind Regards,
>>> James
>>>
>>> Simon Pickett wrote:
>>>
 Hi James,

 What you really need to do is to check out the many freely

>>> available
>>>
 pdfs for R beginners. Here is a good place to start

 http://cran.r-project.org/other-docs.html

 If I am right interpreting what you want, I think you need

>>> to create a
>>>
 blank plot with no axes, axis labels etc. Try

 plot(x,y,xlab="",ylab="",xaxt=NULL,yaxt=NULL,type="n")
 #blank plot
 points(x,y)

 type "?par" into R and see how you can set parameters like

>>> this up as
>>>
 the default.

 Hope this helps?

 Simon.


 - Original Message - From: "James Nicolson"
 
 To: 
 Sent: Sunday, February 15, 2009 10:29 PM
 Subject: [R] Unadulterated plot



> To all,
>
> Apologies if this question has already been asked but I can't find
> anything. I can't seem to think of more specific search
>
>>> terms. I want
>>>
> to display/create a file of a pure plot with a specific height and
> width. I want to utilise every single pixel inside the
>
>>> axes. I do not
>>>
> want to display any margins, legends, axes, titles or
>
>>> spaces around
>>>
> the edges. Is this possible? Additionally, the plot I am
>
>>> working with
>>>
> is a filled.contour plot and I can not remove the legend?
>
>>> How can I
>>>
> do this?
>
> Kind Regards,
> James
>
> __
> 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.
>>>
>>>
>>
>> The contents of this e-mail are confidential and may be subject to legal 
>> privilege.
>>  If you are not the intended recipient you must not use, disseminate, 
>> distribute or
>>  reproduce all or any part of this e-mail or attachments.  If you have 
>> received this
>>  e-mail in error, please notify the sender and delete all material 
>> pertaining to this
>>  e-mail.  Any opinion or views expressed in this e-mail are those of the 
>> individual
>>  sender and may not represent those of The New Zealand Institute for Plant 
>> and
>>  Food Research Limited.
>>
>>
>
> __
> R-help@r-project.o

Re: [R] : record which entry in one file doesn't appear in a different file

2009-02-25 Thread Patrizio Frederic
hey Laura,
I hope this help

f1 = c("a","b","c")
f2 = c("b","a","c","d")
match(f2,f1)
f3 = match(f2,f1,0)
?match
cbind(f2,f3)
cbind(f2,f3>0)
f4 = ifelse(f3>0,"yes","no")
cbind(f2,f4)
data.frame(f2,f4)

Patrizio

2009/2/25 Laura Rodriguez Murillo :
> Hi dear list,
>
> If anybody could help me, it would be great!
>
> I have two files:
> File 1 is a list (one column and around 10 rows)
> File 2 is a list with all the names from file one and a few more (one
> column and more than 10 rows)
>
> What I want is to add a column in file 2 that says which name appeared
> in file 1 and which doesn't (yes and no would work as a code)
> It's very important to keep the order of the names in file 2.
>
> Thank you!
>
> Laura
>
> __
> 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] Way to rotate a histogram?

2009-03-17 Thread Patrizio Frederic
Jason,
be carefully to the order of intensities (counts or densities):

x=rnorm(1000)
par(mfrow=c(2,2))
h=hist(x,breaks=bk<-c(-5,-3,-2,-1,-.5,0,1,3,5))
barplot(rev(h$intensities),rev(bk[2:9]-bk[1:8]),space=0,horiz=T) # compare to
axis(2)
barplot(h$intensities,bk[2:9]-bk[1:8],space=0,horiz=T)
axis(2)

Patrizio

2009/3/17 David Winsemius :
> I believe that hist will return a vector that could be passed to barplot:
>
>  h.islands <- hist(islands)
>
>> barplot(h.islands$intensities, horiz=TRUE)  # or
>> barplot(h.islands$counts, horiz=TRUE)
>
> David Winsemius
>
> On Mar 17, 2009, at 11:38 AM, Jason Rupert wrote:
>
>>
>>
>> Here is what I have so far:
>>>
>>> test_data<-rnorm(100)
>>> par(mfrow=c(1,3)) # 3 rows by 1 columns layout of plots
>>> hist(test_data)
>>> boxplot(test_data)
>>> qqnorm(test_data)
>>
>> I noticed that I can rotate a boxplot via "horizontal", but apparently
>> "hist" does not have that functionality.
>>
>> I tried stacking the plots vertically:
>> test_data<-rnorm(100)
>> par(mfrow=c(3,1)) # 3 rows by 1 columns layout of plots
>> hist(test_data)
>> boxplot(test_data, horizontal=TRUE)
>> qqnorm(test_data)
>>
>> However, I would have to rotate the QQnorm plot, which would be pretty
>> confusing and I think non-standard.
>>
>> Thank you again for any feedback and insight regarding trying to reproduce
>> the JMP figure shown at:
>> http://n2.nabble.com/Can-R-produce-this-plot--td2489288.html
>>
>> --- On Tue, 3/17/09, Jason Rupert  wrote:
>>
>>> From: Jason Rupert 
>>> Subject: Re: [R] R package to automatically produce combination plot?
>>> To: R-help@r-project.org
>>> Date: Tuesday, March 17, 2009, 9:39 AM
>>> I guess no reply means there is not an existing package to
>>> produce the plot?
>>>
>>> I will post the results of my script to hopefully help
>>> others who are trying to formulate the same plot.
>>>
>>> Thanks again.
>>>
>>>
>>> --- On Mon, 3/16/09, Jason Rupert
>>>  wrote:
>>>
 From: Jason Rupert 
 Subject: [R] R package to automatically produce
>>>
>>> combination plot?

 To: R-help@r-project.org
 Date: Monday, March 16, 2009, 8:14 PM
 By any chance is there an R package that automatically
 produces the plot shown at the following link:

>>> http://n2.nabble.com/Can-R-produce-this-plot--td2489288.html

 That is an R package to produce on plot that has the
 following:
 (a) a vertically oriented histogram,
 (b) associated barplot, and
 (c) quantile-quantile plot (Q-Q Plot).

 This is based on a class lecture from University of
 Pennsylvania:
 stat.wharton.upenn.edu/~mcjon/stat-431/lecture-02.pdf

 I am pretty confident I can put one together, but just
 wanted to check that there does not already exist an R
 package to output such a plot.

 Thanks again.

 __
 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-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.
>
> David Winsemius, MD
> Heritage Laboratories
> West Hartford, CT
>
> __
> 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] multiple barplot

2009-03-18 Thread Patrizio Frederic
Hi Agus,
try this two ones

d<-matrix(rpois(45,3),5,9)
barplot(d,beside=T,col=rainbow(5),names=c("CRTL","LSB","ONEMKR",
"TWOMKR","BLUP","BLUPQ","BLUP1M","BLUP2M","GAS"),las=2)


barplot(d,beside=T,col=rainbow(5),names=c("CRTL","LSB","ONEMKR",
"TWOMKR","BLUP","BLUPQ","BLUP1M","BLUP2M","GAS"),cex.names=.8)

best

Patrizio

2009/3/18 Agus Susanto :
> Dear all,
> I want to put 9 barplots side by side. My code below only print 5 names from
> 9 names I gave.
> Problem: how to print all of those 9 names? I use cex=0.8 but did not work,
> it gave me error message.
>
> d<-matrix(rpois(45,3),5,9)
> barplot(d,beside=T,col=rainbow(5),names=c("CRTL","LSB","ONEMKR",
>        "TWOMKR","BLUP","BLUPQ","BLUP1M","BLUP2M","GAS"))
>
> # with cex=0.8 (but not working)
> barplot(d,beside=T,col=rainbow(5),names=c("CRTL","LSB","ONEMKR",
>        "TWOMKR","BLUP","BLUPQ","BLUP1M","BLUP2M","GAS"),cex=0.8)
>
> Thanks in advance.
> --
> Agus Susanto
>
>        [[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-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] repeated searching of no-missing values

2008-12-10 Thread Patrizio Frederic
hi all,
I have a data frame such as:

1 blue  0.3
1 NA0.4
1 red   NA
2 blue  NA
2 green NA
2 blue  NA
3 red   0.5
3 blue  NA
3 NA1.1

I wish to find the last non-missing value in every 3ple: ie I want a 3
by 3 data.frame such as:

1 red   0.4
2 blue  NA
3 blue  1.1

I have written a little script

data = structure(list(V1 = c(1L, 1L, 1L, 2L, 2L, 2L, 3L, 3L, 3L
), V2 = structure(c(1L, NA, 3L, 1L, 2L, 1L, 3L, 1L, NA), .Label = c("blue",
"green", "red"), class = "factor"), V3 = c(0.3, 0.4, NA, NA,
NA, NA, 0.5, NA, 1.1)), .Names = c("V1", "V2", "V3"), class =
"data.frame", row.names = c(NA,
-9L))

cl= function(x) x[max(which(!is.na(x)))]
choose.last = function(x) tapply(x,x[,1],cl)

# now function choose.last works properly on numeric vectors:

> choose.last(data[,3])
  1   2   3
0.4  NA 1.1

# but not on factors (I loose the factor labels):

> choose.last(data[,2])
1 2 3
3 1 1

# moreover, if I apply this function to the whole data.frame
# the output is a character matrix

> apply(data,2,choose.last)
  V1  V2 V3
1 "1" "red"  "0.4"
2 "2" "blue" NA
3 "3" "blue" "1.1"

# and if I sapply, I loose factors labels

> sapply(data,choose.last)
  V1 V2  V3
1  1  3 0.4
2  2  1  NA
3  3  1 1.1

any hint?

Thanks in advance,

Patrizio

+-
| Patrizio Frederic, PhD
| Research associate in Statistics,
| Department of Economics,
| University of Modena and Reggio Emilia,
| Via Berengario 51,
| 41100 Modena, Italy
|
| tel:  +39 059 205 6727
| fax:  +39 059 205 6947
| mail: [EMAIL PROTECTED]
+-

__
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] repeated searching of no-missing values

2008-12-11 Thread Patrizio Frederic
dear Hadley and Bert,
thank you very much for your suggestions. I asked one question and I
learned 2 things:

1. Hadley,

library(plyr)
ddply(data, .(V1), colwise(cl))

that is exactly what I was searching for.

2. Bert,

> ?tapply says that the first argument is an **atomic** vector. A
> factor is not an atomic vector. So tapply interprets it as such by looking
> only at its representation, which is as integer values.

I admit I payed not enough attention to the definition of **atomic** vector.
That implies a deeper understanding of structures of data. I'm working with!

Many thanks,

Patrizio

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