Re: [R] Sensitivity and Specificity

2022-10-24 Thread Prof. Dr. Matthias Kohl

MKclass::perfMeasures(predict_testing, truth = testing$case, namePos = 1)

should also work and computes 80 performance measures.

Best Matthias

Am 25.10.22 um 06:42 schrieb Jin Li:

Hi Greg,

This can be done by:
spm::pred.acc(testing$case,  predict_testing)

It will return both sensitivity and specificity, along with a few other
commonly used measures.

Hope this helps,
Jin

On Tue, Oct 25, 2022 at 6:01 AM Rui Barradas  wrote:


Às 16:50 de 24/10/2022, greg holly escreveu:

Hi Michael,

I appreciate your writing. Here are what I have after;


predict_testing <- ifelse(predict > 0.5,1,0)

head(predict)

   1  2  3  5  7  8
0.29006984 0.28370507 0.10761993 0.02204224 0.12873872 0.08127920


# Sensitivity and Specificity





sensitivity<-(predict_testing[2,2]/(predict_testing[2,2]+predict_testing[2,1]))*100

Error in predict_testing[2, 2] : incorrect number of dimensions

sensitivity

function (data, ...)
{
  UseMethod("sensitivity")
}








specificity<-(predict_testing[1,1]/(predict_testing[1,1]+predict_testing[1,2]))*100

Error in predict_testing[1, 1] : incorrect number of dimensions

specificity

function (data, ...)
{
  UseMethod("specificity")
}



On Mon, Oct 24, 2022 at 10:45 AM Michael Dewey 
wrote:


Rather hard to know without seeing what output you expected and what
error message you got if any but did you mean to summarise your variable
predict before doing anything with it?

Michael

On 24/10/2022 16:17, greg holly wrote:

Hi all R-Help ,

After partitioning my data to testing and training (please see below),

I

need to estimate the Sensitivity and Specificity. I failed. It would be
appropriate to get your help.

Best regards,
Greg


inTrain <- createDataPartition(y=data$case,
  p=0.7,
  list=FALSE)
training <- data[ inTrain,]
testing  <- data[-inTrain,]

attach(training)
#model training and prediction
data_training <- glm(case ~ age+BMI+Calcium+Albumin+meno_1, data =
training, family = binomial(link="logit"))

predict <- predict(data_training, data_predict = testing, type =

"response")


predict_testing <- ifelse(predict > 0.5,1,0)

# Sensitivity and Specificity





  
sensitivity<-(predict_testing[2,2]/(predict_testing[2,2]+predict_testing[2,1]))*100

sensitivity





  
specificity<-(predict_testing[1,1]/(predict_testing[1,1]+predict_testing[1,2]))*100

specificity

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



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



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


Hello,

Instead of computing by hand, why not use package caret?


tbl <- table(predict_testing, testing$case)
caret::sensitivity(tbl)
caret::specificity(tbl)


Hope this helps,

Rui Barradas

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






--
Prof. Dr. Matthias Kohl
www.stamats.de

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


Re: [R] confidence intervals for the difference between group means

2020-08-04 Thread Prof. Dr. Matthias Kohl

you could try:

library(MKinfer)
meanDiffCI(a, b, boot = TRUE)

Best
Matthias

Am 04.08.20 um 16:08 schrieb varin sacha via R-help:

Dear R-experts,

Using the bootES package I can easily calculate the bootstrap confidence 
intervals of the means like in the toy example here below. Now, I am looking 
for the confidence intervals for the difference between group means. In my 
case, the point estimate of the mean difference is 64.4. I am looking at the 
95% confidence intervals around this point estimate (64.4).

Many thanks for your response.


library(bootES)
a<-c(523,435,478,567,654)
b<-c(423,523,421,467,501)
bootES(a)
bootES(b)


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



--
Prof. Dr. Matthias Kohl
www.stamats.de

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


Re: [R] ncol() vs. length() on data.frames

2020-03-31 Thread Prof. Dr. Matthias Kohl

should have added: dim(x)[2L] -> length(x)

Am 31.03.20 um 16:21 schrieb Prof. Dr. Matthias Kohl:

Dear Ivan,

if I enter ncol in the console, I get

function (x)
dim(x)[2L]



indicating that function dim is called. Function dim has a method for 
data.frame; see methods("dim").


The dim-method for data.frame is

dim.data.frame
function (x)
c(.row_names_info(x, 2L), length(x))



Hence, it calls length on the provided data.frame. In addition, some 
"magic" with .row_names_info is performed, where


base:::.row_names_info
function (x, type = 1L)
.Internal(shortRowNames(x, type))



Best
Matthias

Am 31.03.20 um 16:10 schrieb Ivan Calandra:

Thanks Ivan for the answer.

So it confirms my first thought that these two functions are equivalent
when applied to a "simple" data.frame.

The reason I was asking is because I have gotten used to use length() in
my scripts. It works perfectly and I understand it easily. But to be
honest, ncol() is more intuitive to most users (especially the novice)
so I was thinking about switching to using this function instead (all my
data.frames are created from read.csv() or similar functions so there
should not be any issue). But before doing that, I want to be sure that
it is not going to create unexpected results.

Thank you,
Ivan

--
Dr. Ivan Calandra
TraCEr, laboratory for Traceology and Controlled Experiments
MONREPOS Archaeological Research Centre and
Museum for Human Behavioural Evolution
Schloss Monrepos
56567 Neuwied, Germany
+49 (0) 2631 9772-243
https://www.researchgate.net/profile/Ivan_Calandra

On 31/03/2020 16:00, Ivan Krylov wrote:

On Tue, 31 Mar 2020 14:47:54 +0200
Ivan Calandra  wrote:


On a simple data.frame (i.e. each element is a vector), ncol() and
length() will give the same result.
Are they just equivalent on such objects, or are they differences in
some cases?

I am not aware of any exceptions to ncol(dataframe)==length(dataframe)
(in fact, ncol(x) is dim(x)[2L] and ?dim says that dim(dataframe)
returns c(length(attr(dataframe, 'row.names')), length(dataframe))), but
watch out for AsIs columns which can have columns of their own:

x <- data.frame(I(volcano))
dim(x)
# [1] 87  1
length(x)
# [1] 1
dim(x[,1])
# [1] 87 61




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





--
Prof. Dr. Matthias Kohl
www.stamats.de

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


Re: [R] ncol() vs. length() on data.frames

2020-03-31 Thread Prof. Dr. Matthias Kohl

Dear Ivan,

if I enter ncol in the console, I get

function (x)
dim(x)[2L]



indicating that function dim is called. Function dim has a method for 
data.frame; see methods("dim").


The dim-method for data.frame is

dim.data.frame
function (x)
c(.row_names_info(x, 2L), length(x))



Hence, it calls length on the provided data.frame. In addition, some 
"magic" with .row_names_info is performed, where


base:::.row_names_info
function (x, type = 1L)
.Internal(shortRowNames(x, type))



Best
Matthias

Am 31.03.20 um 16:10 schrieb Ivan Calandra:

Thanks Ivan for the answer.

So it confirms my first thought that these two functions are equivalent
when applied to a "simple" data.frame.

The reason I was asking is because I have gotten used to use length() in
my scripts. It works perfectly and I understand it easily. But to be
honest, ncol() is more intuitive to most users (especially the novice)
so I was thinking about switching to using this function instead (all my
data.frames are created from read.csv() or similar functions so there
should not be any issue). But before doing that, I want to be sure that
it is not going to create unexpected results.

Thank you,
Ivan

--
Dr. Ivan Calandra
TraCEr, laboratory for Traceology and Controlled Experiments
MONREPOS Archaeological Research Centre and
Museum for Human Behavioural Evolution
Schloss Monrepos
56567 Neuwied, Germany
+49 (0) 2631 9772-243
https://www.researchgate.net/profile/Ivan_Calandra

On 31/03/2020 16:00, Ivan Krylov wrote:

On Tue, 31 Mar 2020 14:47:54 +0200
Ivan Calandra  wrote:


On a simple data.frame (i.e. each element is a vector), ncol() and
length() will give the same result.
Are they just equivalent on such objects, or are they differences in
some cases?

I am not aware of any exceptions to ncol(dataframe)==length(dataframe)
(in fact, ncol(x) is dim(x)[2L] and ?dim says that dim(dataframe)
returns c(length(attr(dataframe, 'row.names')), length(dataframe))), but
watch out for AsIs columns which can have columns of their own:

x <- data.frame(I(volcano))
dim(x)
# [1] 87  1
length(x)
# [1] 1
dim(x[,1])
# [1] 87 61




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



--
Prof. Dr. Matthias Kohl
www.stamats.de

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


Re: [R] How to rum Multiple ANOVA and Multiple T-test between the same groups?

2019-02-10 Thread Prof. Dr. Matthias Kohl
Have a look at Bioconductor package genefilter, especially functions 
colttests and colFtests.

Best Matthias

Am 10.02.19 um 10:35 schrieb AbouEl-Makarim Aboueissa:

Dear All: good morning





*Re:* How to rum Multiple ANOVA and Multiple T-test between the same groups.



Your help will be highly appreciated.





*1.*  is there a way to run multiple t-tests on different variables between
the same two groups.





*Data for t-tests:*



The data frame “dataTtest”  has 5 variables (x1,x2,x3,x4,x5) and one factor
(factor1) with 2 levels (group1, group2).





x1<-rnorm(20,1,1)

x2<-rnorm(20,2,1)

x3<-rnorm(20,3,1)

x4<-rnorm(20,4,1)

x5<-rnorm(20,5,1)

factor1<-rep(c("group1", "group2"), each = 10)

dataTtest<-data.frame(x1,x2,x3,x4,x5,factor1)

dataTtest









*2.* is there a way to run *multiple ANOVA* and multiple comparisons *Tukey
tests* on different variables between the same groups.





*Data for ANOVA tests:*



The data frame “dataANOVA”  has 6 variables (x1,x2,x3,x4,x5,x6) and one
factor (factor2) with 5 levels (group1, group2, group3, group4, group5).







x1<-rnorm(40,1,1)

x2<-rnorm(40,2,1)

x3<-rnorm(40,3,1)

x4<-rnorm(40,4,1)

x5<-rnorm(40,5,1)

x6<-rnorm(40,6,1)

factor2<-rep(c("group1", "group2", "group3", "group4", "group5"), each = 8)

dataANOVA<-data.frame(x1,x2,x3,x4,x5,x6,factor2)

dataANOVA





with many thanks

abou
__


*AbouEl-Makarim Aboueissa, PhD*

*Professor, Statistics and Data Science*
*Graduate Coordinator*

*Department of Mathematics and Statistics*
*University of Southern Maine*

[[alternative HTML version deleted]]

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



--
Prof. Dr. Matthias Kohl
www.stamats.de

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


Re: [R] Paired sample t-test with mi.t.test

2017-04-07 Thread Prof. Dr. Matthias Kohl

Dear Joel,

are you trying to apply function mi.t.test from my package MKmisc?

Could you please try:
mi.t.test(implist, "pre_test", "post_test", alternative =
"greater", paired = TRUE, var.equal = TRUE, conf.level = 0.95)

x and y are the names of the variables, not the variables themselves.

Best
Matthias

Am 06.04.2017 um 18:32 schrieb Joel Gagnon:

Dear all,

It is my first time posting on this list so forgive me for any rookie
mistakes I could make.

I want to conduct t-tests on a dataset that has been imputed using the mice
package:
imput_pps <- mice(pps, m=20, maxit=20, meth='pmm') # pps is my dataset. It
contains items from an 11-item questionnaire gather at pre and post test.
So the data set has 22 columns.

I then proceed to compute the total scores for the pre and post test on my
imputed datasets:

long_pps <- complete(imput_pps, action ="long", include = TRUE)
long_pps$pre_test <- rowSums(long_pps[ ,c(3:13)])
long_pps$post_test <- rowSums(long_pps[ , c(14:24)])

I then used as.mids to convert back to mids object:
mids_pps <- as.mids(long_pps)

Next, I created an imputation list object using mitools:
implist <- lapply(seq(mids_pps$m), function(im) complete(mids_pps, im))
implist <- imputationList(implist)

Now, I want to conduct t-tests using the mi.t.test package. I tried the
following code:
mi.t.test(implist, implist$pre_test, implist$post_test, alternative =
"greater", paired = TRUE, var.equal = TRUE, conf.level = 0.95)

When I run this code, R tells me that Y is missing. I know this may sound
stupid, but I thought that I specified Y with this line: implist$pre_test,
implist$post_test - with implist$pre_test being X and implist$post_test
being Y - like I usually do for a normal t-test using the t.test function.

It seems I don't quite understand what the Y variable is supposed to
represent. Could someone help me figure out what I am doing wrong? You
help would be very much appreciated.

Best regards,

Joel Gagnon, Ph.D(c),
Department of Psychology,
Université du Québec à Trois-Rivières
Québec, Canada

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



--
Prof. Dr. Matthias Kohl
www.stamats.de

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

Re: [R] How to install pkg when "not found" ?

2015-12-30 Thread Prof. Dr. Matthias Kohl

use quotes!
install.packages("ggplot2")

Am 30.12.2015 um 06:41 schrieb Judson:


Using "Install Packages" from CRAN, in RStudio on Windows 7,
I downloaded (and supposedly installed)  ggplot2 package to here:

C:\Program Files\R\R-3.1.0\library\ggplot2_2.0.0\ggplot2


 when I try this:

require(ggplot2)

 I get the following:
Loading required package: ggplot2
Warning message:
In library(package, lib.loc = lib.loc, character.only = TRUE, logical.return = 
TRUE,  :

   there is no package called 'ggplot2'


 I also tried this:

install.packages(ggplot2)

Error in install.packages : object 'ggplot2' not found

..


The above library has the usual things like MASS, stats, base, utils, graphics 
and so on
. and if I do this:

require(MASS)

 I get this:
Loading required package: MASS

 so that works


. Any idea what I'm doing wrong with ggplot2?
... judson blake

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



--
Prof. Dr. Matthias Kohl
www.stamats.de

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


Re: [R] C.D.F

2014-08-11 Thread Prof. Dr. Matthias Kohl

Dear Diba,

you could try package distr; eg.

library(distr)
G1 - Gammad(scale = 0.7, shape = 0.5)
G2 - Gammad(scale = 2.1, shape = 1.7)
G3 - G1+G2 # convolution
G3

For the convolution exact formulas are applied if available, otherwise 
we use FFT; see also http://www.jstatsoft.org/v59/i04/ (will appear 
soon) resp. a previous version at http://arxiv.org/abs/1006.0764


hth
Matthias

Am 11.08.2014 um 10:17 schrieb pari hesabi:

Hello everybody,

Can anybody help me to write a program for the CDF of sum of two independent 
gamma random  variables ( covolution of two gamma distributions)  with 
different amounts of parameters( the shape parameters are the same)?

Thank you

Diba

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



--
Prof. Dr. Matthias Kohl
www.stamats.de

__
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] legendre quadrature

2014-05-01 Thread Prof. Dr. Matthias Kohl

you could use package distrEx:

library(distrEx)
GLIntegrate(function(x) x^2, lower = -1, upper = 1, order = 50)

hth
Matthias

On 01.05.2014 09:43, pari hesabi wrote:

Hello everybody
I need to approximate the amount of integral by using
  legendre quadrature. I have written a program which doesn't give me a
logical answer; Can anybody help me and send the correct program? For
example  the approximated amount of integral of ( x ^2)  on (-1,1) based
  on legendre quad rule.




integrand-function(x) {x^2}
rules - legendre.quadrature.rules( 50 )
order.rule - rules[[50]]
chebyshev.c.quadrature(integrand, order.rule, lower = -1, upper = 1)

Thank you
Diba
[[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.



--
Prof. Dr. Matthias Kohl
www.stamats.de

__
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] joint distribution

2014-04-22 Thread Prof. Dr. Matthias Kohl

have a look at our package distr:

library(distr)
x1 - Norm(mean = 0, sd = 1)
x2 - Binom(size = 1, prob = 0.75)
x3 - x1 + x2
plot(x3)
# to get density, cdf, quantile, and random numbers use
d(x3)(5)
p(x3)(0)
q(x3)(0.975)
r(x3)(20)

# you can also have additonal coefficients; eg.
x4 - 3*x1 + 2*x2
plot(x4)

For more details on the computations see http://arxiv.org/abs/1006.0764

hth,
Matthias

On 22.04.2014 13:11, Ms khulood aljehani wrote:


Hello
i have two independent variablesx1 from normal (0,1)x2 from bernoulli (o.75)
i need the density estimation of(b1*x1) + (b2*x2)
where b1 and b2 are two fixed coefficients
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.



--
Prof. Dr. Matthias Kohl
www.stamats.de

__
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] Creating new instances from original ones

2013-03-31 Thread Prof. Dr. Matthias Kohl

see function SMOTE in package DMwR
hth
Matthias

On 31.03.2013 10:46, Nicolás Sánchez wrote:

I have a question about data mining. I have a dataset of 70 instances with
14 features that belong to 4 classes. As the number of each class is not
enough to obtain a good accuracy using some classifiers( svm, rna, knn) I
need to oversampling the number of instances of each class.

I have heard that there is a method to do this. It consists in generating
these new instances as follows:

new_instance  original_instance + u(epsilon)

  U(epsilon) is a uniform number in the range [-epsilon,epsilon] and this
number is applied to each feature of the dataset to obtain a new instance
without modified the original class.

Anybody has used this method to oversampling his data? Anybody has more
information about it?

Thanks in advance!

[[alternative HTML version deleted]]

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



--
Prof. Dr. Matthias Kohl
www.stamats.de

__
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] integrate function

2013-02-12 Thread Prof. Dr. Matthias Kohl

use pmin instead of min.
hth
Matthias

On 12.02.2013 16:41, dan wang wrote:

Hi All,

Can any one help to explain why min and max function couldn't work in the
integrate function directly.

For example, if issue following into R:

integrand - function(x) {min(1-x, x^2)}
integrate(integrand, lower = 0, upper = 1)

it will return this:
Error in integrate(integrand, lower = 0, upper = 1) :
   evaluation of function gave a result of wrong length


However, as min(U,V) = (U+V)/2-abs(U-V)/2

Below code working;

integrand - function(x){(1-x+x^2)/2-0.5*abs(1-x-x^2)}

integrate(integrand, lower = 0, upper = 1)

0.151639 with absolute error  0.00011

I am confused...


Dan

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



--
Prof. Dr. Matthias Kohl
www.stamats.de

__
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] Convolve 2 uniform distributed variables

2012-05-24 Thread Prof. Dr. Matthias Kohl

take a look at the distr package

library(distr)
U - Unif()
U2 - U + U
# or
U2 - convpow(U, 2)
plot(U2)
# or
curve(d(U2)(x), from = 0, to = 2)


Best
Matthias

On 24.05.2012 08:29, c...@mail.tu-berlin.de wrote:

Hi,
I want to convolve 2 uniform distributed [0,1] variables, so that I get
this graphic:
http://de.wikipedia.org/wiki/Benutzer:Alfred_Heiligenbrunner/Liste_von_Verteilungsdichten_der_Summe_gleichverteilter_Zufallsvariabler
(second graphic)

if I do

u1-seq(0,1,length=100)
fu1=dunif(u1,min=0,max=1)

plot(u1,fu1,type=l,xlim=c(-2,2))


u2-seq(0,1,length=100)
fu2=dunif(u2,min=0,max=1)

u1u2-convolve(u1,rev(u1),typ=o)

plot(u1u2)

it does not work, could you help me please? The point is the convolve
function I think, what do I have to type, to get the correct convolution
of two uniform distributed [0,1] variables as to be seen in the second
graphic in the given link?

Thanks a lot

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


--
Prof. Dr. Matthias Kohl
www.stamats.de

__
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] Compare String Similarity

2012-04-19 Thread Prof. Dr. Matthias Kohl

you should also look at Bioconductor Package Biostrings
hth
Matthias

Am 19.04.2012 18:01, schrieb R. Michael Weylandt:

Though if you do decide to use Levenstein, it's implemented here in R:
http://finzi.psych.upenn.edu/R/library/RecordLinkage/html/strcmp.html

I'm pretty sure this is in C code so it should be mighty fast.

Michael

On Thu, Apr 19, 2012 at 11:40 AM, Bert Guntergunter.ber...@gene.com  wrote:

Wrong list.This is R, not statistics (or linguistics?).Please post elsewhere.

-- Bert

On Thu, Apr 19, 2012 at 8:05 AM, Alekseiy Beloshitskiy
abeloshits...@velti.com  wrote:

Dear All,

I need to estimate the level of similarity of two strings. For example:
string1- c(depending,audience,research, school);
string2- c(audience,push,drama,button,depending);

The words in string may occur in different order though. What function would 
you recommend to use to estimate similarity (e.g., levenstein, distance)?

Appreciate for any advices.

-Alex

[[alternative HTML version deleted]]

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



--

Bert Gunter
Genentech Nonclinical Biostatistics

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

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

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


--
Prof. Dr. Matthias Kohl
www.stamats.de

__
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] Sweave UFT8 problem

2012-04-15 Thread Prof. Dr. Matthias Kohl

try:
Sweave(sim_pi.Rnw, encoding = utf8)

Best,
Matthias

On 15.04.2012 11:41, Philippe Grosjean wrote:

Hello,

Have you tried to put that command in a comment:

%\usepackage[utf8]{inputenc}

I haven't tested it in this particular case, but it works in some other
situations.
Best,

Philippe

..¡}))
) ) ) ) )
( ( ( ( ( Prof. Philippe Grosjean
) ) ) ) )
( ( ( ( ( Numerical Ecology of Aquatic Systems
) ) ) ) ) Mons University, Belgium
( ( ( ( (
..

On 14/04/12 22:37, Mark Heckmann wrote:

Hi,

I work on MacOS, trying to Sweave an UFT8 document.
AFAI remember R 2.14 used to render a warning when the encoding was not
declared when using Sweave.
With R 2.15 it seems to render an error.

Sweave(sim_pi.Rnw)
Error: 'sim_pi.Rnw' is not ASCII and does not declare an encoding

Declaring an encoding by adding a line like

\usepackage[utf8]{inputenc}

in the preamble does the job. In my case though the .Rnw document does no
have a preamble as it is just one chapter.
All chapters are Sweaved separately (due to computation time). Hence I
cannot inject the above line as LaTex will cause an error afterwards.
(usepackage{} is only allowed in the preamble which only appears once in
the main document, not in each chapter).

How can I get around this not using the terminal for Sweaving, like e.g.
R CMD Sweave --encoding=utf-8 sim_pi.Rnw ?

Thanks
Mark

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


--
Prof. Dr. Matthias Kohl
www.stamats.de

__
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] calculate quantiles of a custom function

2012-01-03 Thread Prof. Dr. Matthias Kohl

Dear Gerhard,

you could also use package distr; e.g.

library(distr)

## use generating function AbscontDistribution
D - AbscontDistribution(d = function(x) dbeta(x, 2, 6) + dbeta(x,6,2), 
low = 0, up = 1, withStand = TRUE)


## quantiles
q(D)(seq(0,1,0.1))

Best
Matthias

On 03.01.2012 19:33, Albyn Jones wrote:

What do quantiles mean here? If you have a mixture density, say

myf - function(x,p0) p0*dbeta(x,2,6) + (1-p0)*dbeta(x,6,2)

then I know what quantiles mean. To find the Pth quantile use uniroot to
solve for the x such that myf(x,p0) - P =0.

albyn

Quoting VictorDelgado victor.m...@fjp.mg.gov.br:



Gerhard wrote



Suppose I create a custom function, consisting of two
beta-distributions:

myfunction - function(x) {
dbeta(x,2,6) + dbeta(x,6,2)
}

How can I calculate the quantiles of myfunction?

Thank you in advance,

Gerhard




Gehard, if do you want to know the quantiles of the new distribution
created
by myfunction. Maybe you can also do:

x - seq(0,1,.01) # insert your 'x'
q - myfunction(x)
# And:
quantile(x)

0% 25% 50% 75% 100%
0.00 1.476177 2.045389 2.581226 2.817425

# This gives the sample quantiles. You can also look foward to
simulations
(like Bert Gunter had suggested) to know better the properties of
distributions quantiles obtained after 'myfunction'.



-
Victor Delgado
cedeplar.ufmg.br P.H.D. student
www.fjp.mg.gov.br reseacher
--
View this message in context:
http://r.789695.n4.nabble.com/calculate-quantiles-of-a-custom-function-tp4256887p4257551.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.


--
Prof. Dr. Matthias Kohl
www.stamats.de

__
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] Defining functions inside loops

2011-02-15 Thread Dr. Matthias Kohl
.





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


--
Dr. Matthias Kohl
www.stamats.de

__
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] Package distr and define your own distribution

2011-02-12 Thread Dr. Matthias Kohl

Dear Gabriel,

it is possible. You can define a new class or just an object of an 
already existing class.


Did you look at:
library(distr)
vignette(newDistributions)

as well as
?AbscontDistribution
?DiscreteDistribution

Please let me know if you have further questions!

Best
Matthias

On 11.02.2011 16:05, gabriel.ca...@ubs.com wrote:


Hi all
I  am using the Package distr (and related)

Do you know if it is possible to define your own distribution (object)
GIVEN that you have an analytical form of the probability density
function (pdf) ?

I would then like to use the standard feature of the distr and related
packages.


Best regards

Giuseppe Gabriel Cardi




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


--
Dr. Matthias Kohl
www.stamats.de

__
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 generate a random data from a empirical distribition

2010-07-27 Thread Dr. Matthias Kohl

Hi Dennis,

you should take a look at the CRAN task view for distributions
http://cran.r-project.org/web/views/Distributions.html

Beside that our distr-family of packages might be useful, see also
http://www.jstatsoft.org/v35/i10/
http://cran.r-project.org/web/packages/distrDoc/vignettes/distr.pdf

Best,
Matthias

On 27.07.2010 10:37, Dennis Murphy wrote:

Hi:

On Mon, Jul 26, 2010 at 11:36 AM, xin weixin...@stat.psu.edu  wrote:



hi, this is more a statistical question than a R question. but I do want to
know how to implement this in R.
I have 10,000 data points. Is there any way to generate a empirical
probablity distribution from it (the problem is that I do not know what
exactly this distribution follows, normal, beta?). My ultimate goal is to
generate addition 20,000 data point from this empirical distribution
created
from the existing 10,000 data points.
thank you all in advance.



The problem, it seems to me, is the leap of faith you're taking that the
empirical distribution of your manifest sample will serve as a useful
data-generating mechanism for the 20,000 future observations you want to
take. I would think that, if you intend to take a sample of 20,000 from ANY
distribution, you would want some confidence in the specification of said
distribution.

Even if you don't know exactly what type of population distribution you're
dealing with, there are ways to narrow down the set of possibilities. What
is the domain/support of the distribution? For example, the Normal is
defined on all of R (as in the real numbers, not our favorite statistical
programming language), whereas the lognormal, Gamma and Weibull
distributions, among others, are defined on the nonnegative reals. The beta
distribution is defined on [0, 1]. Therefore, knowledge of the domain is
useful in and of itself. Is it plausible that the distribution is symmetric,
or should it have a distinct left or right skew? (Similar comments apply to
discrete distributions.) Is censoring or truncation a relevant concern? If
there is a random process that well describes how the data you observe are
generated, that will certainly narrow down the class of potential
data-generating mechanisms/distributions.

Once you've narrowed down the class of possible distributions as much as
possible, you could look into the fitdistr() function in MASS or the
fitdistrplus package on CRAN to test out which candidates seem plausible wrt
your existing sample and which are not. You are not likely to be able to
narrow it down to one family of distributions, but you should have a much
better idea about the characteristics of the distribution that gave rise to
your sample of 10,000 (assuming, of course, that it is a *random* sample)
after going through this exercise, which you can apply to the generation of
the next 20,000 observations.

OTOH, if your existing 10,000 observations were not produced by some random
process, all bets are off.

HTH,
Dennis




--
View this message in context:
http://r.789695.n4.nabble.com/how-to-generate-a-random-data-from-a-empirical-distribition-tp2302716p2302716.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.



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


--
Dr. Matthias Kohl
www.stamats.de

__
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] Robust Asymptotic Statistics (RobASt)

2010-06-07 Thread Dr. Matthias Kohl

Dear Alex,

after installation of ROptRegTS there should be a folder scripts in 
the package directory of ROptRegTS which includes some further examples.
If you are assuming normal errors you should switch to package RobRex 
which is optimized for normal errors. There you will also find a 
folder scripts in the package directory.


In general, multivariate regression should work. But, I have to admit 
that we didn't try many multivariate examples so far. Beside that, the 
user interface is still rudimentary at the moment and some ideas which 
will clearly reduce computation time (as in the case of package RobLox) 
are unfortunately not yet implemented.


So, if you are willing to accept these drawbacks I will be glad to give 
you a hand to get our estimators to work for your data.


Best,
Matthias


On 06.06.2010 06:47, Alex Weslowski wrote:

Hi all,

Other than:
http://www.stamats.de/F2006.r

Are there other good simple examples out there of using the ROptRegTS
package (part of the RobASt project)? I'm hoping to plug it in for
multivariate regression. Or is this not a good idea? Just trying to
find out how it compares to rlm, lts, glm, etc. Hopefully this makes
sense, I'm new to the world of statistics and R.

Thanks!

St0x

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


--
Dr. Matthias Kohl
www.stamats.de

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


[R] [R-pkgs] new package RFLPtools

2010-02-20 Thread Dr. Matthias Kohl

The new package RFLPtools is available on CRAN.

RFLPtools provides analysis functions for DNA fragment molecular weights 
(e.g.\ derived from RFLP-analysis) and nucleotide sequence similarities. 
It aims mainly at the identification of similar or identical fragment 
patterns to evaluate the amount of different genotypes gained from
environmental samples during diversity studies and at further analysis 
of similarities of nucleotide sequences derived from pairwise sequence 
alignments (e.g.\ derived from standalone BLAST). Additionally, 
functions to check the data quality of molecular fingerprints are 
available. To identify the organisms represented by the extracted 
fragment patterns (e.g.\ RFLP genotypes), RFLPtools includes a function 
to compare samples with fragment patterns stored in a reference dataset, 
from which the taxonomic affiliations are already known. To identify 
unknown samples in scientific projects, DNA sequences are used, and the 
gained DNA sequences are compared to existing sequence data via 
alignments tools, as standalone BLAST.
RFLPtools offers tools to generate groups based on tabular report files 
of sequence comparison of pairwise nucleotide sequence alignment.


To get a first impression try:

vignette(RFLPtools)

as well as

library(RFLPtools)
example(RFLPplot)
example(RFLPrefplot)

Best regards,
Fabienne
Alexandra
Matthias

--
Dr. Matthias Kohl
www.stamats.de

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

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


[R] [R-pkgs] new version of RobASt-family of packages

2009-12-09 Thread Matthias Kohl
The new version 0.7 of our RobASt-family of packages is available on 
CRAN for several days. As there were many changes we will only sketch 
the most important ones here. For more details see the corresponding 
NEWS files (e.g. news(package = RobAStBase) or using function NEWS 
from package startupmsg i.e. NEWS(RobAStBase)).



First of all the new package RobLoxBioC was added to the family which 
includes S4 classes and methods for preprocessing omics data, in 
particular gene expression data.



##
## All packages (RandVar, RobAStBase, RobLox, RobLoxBioC,
## ROptEst, ROptEstOld, ROptRegTS, RobRex)
##
- TOBEDONE file was added as a starting point for collaborations. The
file can be displayed via function TOBEDONE from package startupmsg;
e.g. TOBEDONE(distr)
- tests/Examples folder for some automatic testing was introduced


##
## Package RandVar
##
- mainly fixing of warnings and bugs


##
## Package RobAStBase
##
- enhanced plotting, in particular methods for qqplot
- unified treatment of NAs
- extended implementation for total variation neighbourhoods
- implementation of k-step estimator construction extended


##
## Package RobLox
##
- na.rm argument added
- introduction of finite-sample correction


##
## Package ROptEst
##
- optional use of alternative algorithm to obtain Lagrange multipliers
 using duality based optimization
- extended implementation for total variation neighbourhoods
- solutions for general parameter transformations with nuisance
 components
- several extensions to the examples in folder scripts
- implementation of k-step estimator construction extended


##
## Package ROptEstOld
##
- still needed for packages ROptRegTS and RobRex
- removed Symmetry and DistributionSymmetry implementation to
 make ROptEstOld compatible with distr 2.2


##
## Package ROptRegTS
##
- still depends on ROptEstOld


##
## Packages RobRex
##
- moved some of the examples in \dontrun{} to reduce check time ...
- some minor corrections in ExamplesEstimation.R in folder scripts


Best
Peter
Matthias

--
Dr. Matthias Kohl
www.stamats.de

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

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


[R] [R-pkgs] new version of distr-family of packages

2009-11-24 Thread Matthias Kohl
The new version 2.2 of our distr-family of packages has now been 
available on CRAN for several days. As there were many changes we will 
only sketch the most important ones here. For more details see the 
corresponding NEWS files (e.g. news(package = distr) or using function 
NEWS from package startupmsg i.e. NEWS(distr)).


First of all a new package distrEllipse was added to the family which 
includes S4 classes and methods for elliptical symmetric distributions.


##
## All packages (distr, distrEx, distrSim, distrTEst, distrTeach,
## distrMod, distrEllipse, distrDoc as well as
## startupmsg and SweaveListingUtils)
##
- TOBEDONE file was added as a starting point for collaborations. The
 file can be displayed via function TOBEDONE from package startupmsg;
 e.g. TOBEDONE(distr)
- tests/Examples folder for some automatic testing was introduced


##
## Package distr
##
- new diagnostic function (more specifically S4 method)
 qqplot()
 for quantile-quantile plots to work with distribution objects
 as argument;  in addition to function qqplot() from package stats
 has functionality for point-wise and simultaneous confidence bands
- under the hood: a new slot symmetry was introduced
 for distributions

##
## Package distrEx
##
- enhanced E (expectation), m1df, m2df (clipped 1. and 2.
 moments) methods  --- now you can specify lower and
 upper integration bounds in the E(.), var(.) functions
- new class GPareto for generalized Pareto distribution
 (implemented by Nataliya Horbenko)

##
## Package distrMod
##
- methods for qqplot() for parametric model objects as argument
- unified treatment of NAs

##
## Package SweaveListingUtils
##
- commands:
* SweaveListingPreparations() is more flexible
   [see help  NEWS file]
* integrated Andrew Ellis's nice ideas into SweaveListingUtils
   to use \lstenvironment
* individual markup style for base or recommended packages
  (checked for with new function isBaseOrRecommended()) is
   distinct now by default (extra color RRecomdcolor)
* temporarily changes (like background color) made easier:
   by new functions
+ lstdefRstyle to change Rstyle
+ lstsetRall to change all R-like styles (i.e. Rstyle, Rinstyle,
   Routstyle, Rcodestyle)
* now can specify markup styles for Sinput, Soutput, Scode,
   separately as Rin, Rout, Rcode
* colors now have suffix color, i.e.
Rcomment - Rcommentcolor, Rout - Routcolor
- vignette:
* included an example with escape sequences in vignette
* included an example with framed code in vignette

Best
Peter
Matthias

--
Dr. Matthias Kohl
www.stamats.de

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

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


[R] COURSE: Introduction to Metabolomics and biomarker research

2009-10-05 Thread Matthias Kohl

Introduction to Metabolomics and biomarker research.

The course will take place in Innsbruck, 16th to 17th November 2009, and
will be held in German.

For more details see:
http://www.gdch.de/vas/fortbildung/kurse/fortbildung2009.htm#1175

--
Dr. Matthias Kohl
www.stamats.de

__
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] help with median for each row

2009-08-21 Thread Matthias Kohl
if your matrix has many rows you might want to consider rowMedians from 
Bioconductor package Biobase.

hth,
Matthias

Greg Hirson schrieb:

Edward,

See ?apply

x = matrix(rnorm(100), ncol = 10)
apply(x, 1, median)

Hope this helps,

Greg

Edward Chen wrote:

Hi,

I tried looking through google search on whether there's a way to 
computer
the median for each row of a nxn matrix and return the medians for 
each row

for further computation.
And also if the number of columns in the matrix are even, how could I
specify which median to use?

Thank you very much!

  




--
Dr. Matthias Kohl
www.stamats.de

__
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 change the quantile method in bwplot

2009-07-21 Thread Matthias Kohl

for teaching purposes I wrote a corresponding function; cf.
qbxp.stats (as well as qboxplot ...) in package MKmisc.

hth,
Matthias

Deepayan Sarkar schrieb:

On Tue, Jul 21, 2009 at 7:47 AM, Jun Shenjun.shen...@gmail.com wrote:
  

Uwe,

Thank you for your reply.  I am still not very clear about the meanings of
the arguments in the stats function.   To make it clearer, quantile() uses
type=7 as default method. I believe this is the method bwplot() uses to
calculate the quantiles. I want to use type=6 method for bwplot(). How do I
achieve that? Thanks again.



Maybe this will be clearer: bwplot() uses the boxplot.stats() function
to compute the quantiles used, which in turn uses fivenum(), which
has its own quantile calculation (and does not explicitly use
quantile()). There is no easy way to allow for type=6 etc. here.

bwplot() allows you to replace boxplot.stats() and provide your own
alternative. So what you need to do is:

(1) write a function, say, 'my.boxpot.stats', that takes the same
arguments as boxplot.stats() and returns a similar result, but using
your preferred calculation for the quantiles. There are many ways to
do this.

(2) plug in this function into the bwplot() call; e.g. bwplot(...,
stats = my.boxplot.stats)

-Deepayan


  

Jun

2009/7/21 Uwe Ligges lig...@statistik.tu-dortmund.de



Jun Shen wrote:

  

Hi, everyone,

Since quantile calculation has nine different methods in R, I wonder how I
specify a method when calling the bwplot() in lattice. I couldn't find any
information in the documentation. Thanks.




bwplot() uses the panel function panel.bwplot() which allows to specify a
function that calculates the statistics in its argument stats that defaults
to boxplot.stats(). Hence you can change that function.

Example with some fixed values:

bwplot( ~ 1:10,
   stats = function(x, ...)
   return(list(stats=1:5, n=10, conf=1, 10, out=integer(0)))
)


Uwe Ligges

  

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


--
Dr. Matthias Kohl
www.stamats.de

__
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] trouble using optim for maximalisation of 2-parameter function

2009-07-19 Thread Matthias Kohl

try:

# first argument of llik.nor has to be the parameter
llik.nor-function(theta,x){mu-theta[1];sig-theta[2];-length(x)*log(sqrt(2*pi))-length(x)*log(sig)-(1/(2*sig**2))*sum((x-mu)**2)}

# optim by default does minimization
# setting fnscale = -1 one obtains a maximization problem
optim(c(1,5), llik.nor, x=x, control = list(fnscale = -1))

hth,
Matthias

Anjali Sing schrieb:

Hello, I am having trouble using optim.

I want to maximalise a function to its parameters [kind of like: univariate
maximum likelihood estimation, but i wrote the likelihood function myself
because of data issues ]

When I try to optimize a function for only one parameter there is no
problem:

llik.expo-function(x,lam){(length(x)*log(lam))-(length(x)*log(1-exp(-1*lam*
*cons*)))-lam*sum(x)}

cons-

data-c(.)

expomx-optimize(llik.expo,c(0,100),maximum=TRUE,tol=0.0001,x=data)
expomx


 To optimize to two parameters you can't use optimize, so I tried the
following to test my input:

  
llik.nor-function(x,theta){mu-theta[1];sig-theta[2];-length(x)*log(sqrt(2*pi))-length(x)*log(sig)-(1/(2*sig**2))*sum((x-mu)**2)}
  

x-c(-1,4,6,4,2)
normx-optim(c(1,20),llik.nor)



the output should be close to
parameters: mu=3 and sigma=2.366
[This I calculated by hand to compare with the output]

but in stead of output I get an error:
Error in fn(par, ...) : argument theta is missing, with no default

I don't understand why this happened? I hope someone can help me with this
for I am still a [R]ookie.

Kind regards,
Sonko Lady  [Anjali]

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


--
Dr. Matthias Kohl
www.stamats.de

__
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] implementing Maximum Likelihood with distrMod when only the PDF is known

2009-06-24 Thread Matthias Kohl

Dear Guillaume,

retrospectively, I'm not sure if the decision to have spezial initialize 
methods is the optimal way to go. In distrMod and our other packages on 
robust statistics we don't introduce special initialize methods, but use 
so-called generating functions. This approach has the advantage that the 
default initialize method can be used for programming where I find it 
very useful.


I would say it is the canonical (and recommended) way to first define a 
function as generic (like mu) and then implement methods for it.


Choosing names for which there is already an existing definition for a 
generic function may also not be what you want in general. By defining 
new methods you have to respect the definition of the generic function; 
that is, your method definition has to be with respect to the arguments 
in the given generic function and also dispatching will be on these 
arguments (cf. scale in the distrMod example).


In defining our classes we decided that at least the slot r has to be 
filled (of class function) whereas d, p and q are of class 
OptionalFunction. Hence, there are functions to fill d, p and q for 
given r but, so far not for filling r, p and q given d.


A way to avoid implementing r is given in
http://www.stamats.de/distrModExample1.R

I also do not fill the slots p and q in this example. This avoids the 
simulation of a large random sample and speeds up the computation of the 
MLE.


However, this is rather a quick and dirty solution and it would of 
course be better to have a valid definition for r, d, p and q.


Best,
Matthias


guillaume.martin schrieb:

Dear Mathias,
That's pretty amazing, thanks a lot ! I'll have to look all this 
through because I don't easily understand why each part has to be set 
up, in particular the initialize method part seems crucial and is 
not easy to intuit. From what I get, the actual name we give to a 
parameter (my original mu for example) is important in itself, and 
if we introduce new variable names we have to define a new generic, 
right? The simplest option then is to re-use an existing variable name 
that has the same properties/range, right?


Another general question: my actual pdf is of the same type but not 
the exact same as the skew normal. In particular, I don't have a rule 
for building the slot r (eg the one borrowed from the sn package in 
your example); is it a problem? isn't it sufficient to give slot d, 
and then you have automatic methods implemented to get from d() to r() 
slots etc. is that right?


Thanks a lot for your help and time !

Best,

Guillaume


Matthias Kohl a écrit :

Dear Guillaume,

thanks for your interest in the distrMod package.

Regarding your question I took up your example and put a file under:

http://www.stamats.de/distrModExample.R

Hope that helps ...

Don't hesitate to contact me if you have further questions!

Best,
Matthias

guillaume.martin schrieb:

Dear R users and Dear authors of the distr package and sequels

I am trying to use the (very nice) package distrMod as I want to 
implement maximum likelihood (ML) fit of some univariate data for 
which I have derived a theoretical continuous density (pdf). As it 
is a parametric density, I guess that I should implement myself a 
new distribution of class AbscontDistributions (as stated in the pdf 
on creating new distributions in distr), and then use 
MLEstimator() from the distrMod package. Is that correct or is there 
a simpler way to go? Note that I want to use the distr package 
because it allows me to implement simply the convolution of my 
theoretical pdf with some noise distribution that I want to model in 
the data, this is more difficult with fitdistr or mle.
It proved rather difficult for me to implement the new class 
following all the steps provided in the example, so I am asking if 
someone has an example of code he wrote to implement a parametric 
distribution from its pdf alone and then used it with MLEstimator().


I am sorry for the post is a bit long but it is a complicate 
question to me and I am not at all skillful in the handling of such 
notions as S4 - class, etc.. so I am a bit lost here..


As a simple example, suppose my theoretical pdf is the skew normal 
distribution (available in package sn):


#skew normal pdf (default values = the standard normal N(0,1)

fsn-function(x,mu=0,sd=1,d=0) {u = (x-mu)/sd;  f = 
dnorm(u)*pnorm(d*u); return(f/sd)}


# d = shape parameter (any real), mu = location (any real), sd = 
scale (positive real)


#to see what it looks like try
x-seq(-1,4,length=200);plot(fsn(x,d=3),type=l)

#Now I tried to create the classes SkewNorm and 
SkewNormParameter copying the example for the binomial

##Class:parameters
setClass(SkewNormParameter,
representation=representation(mu=numeric,sd=numeric,d=numeric),
prototype=prototype(mu=0,sd=1,d=0,name=gettext(Parameter of the 
Skew Normal distribution)),

contains=Parameter
)

##Class: distribution (created using the pdf of the skew normal 
defined

Re: [R] implementing Maximum Likelihood with distrMod when only the PDF is known

2009-06-23 Thread Matthias Kohl

Dear Guillaume,

thanks for your interest in the distrMod package.

Regarding your question I took up your example and put a file under:

http://www.stamats.de/distrModExample.R

Hope that helps ...

Don't hesitate to contact me if you have further questions!

Best,
Matthias

guillaume.martin schrieb:

Dear R users and Dear authors of the distr package and sequels

I am trying to use the (very nice) package distrMod as I want to 
implement maximum likelihood (ML) fit of some univariate data for 
which I have derived a theoretical continuous density (pdf). As it is 
a parametric density, I guess that I should implement myself a new 
distribution of class AbscontDistributions (as stated in the pdf on 
creating new distributions in distr), and then use MLEstimator() 
from the distrMod package. Is that correct or is there a simpler way 
to go? Note that I want to use the distr package because it allows me 
to implement simply the convolution of my theoretical pdf with some 
noise distribution that I want to model in the data, this is more 
difficult with fitdistr or mle.
It proved rather difficult for me to implement the new class following 
all the steps provided in the example, so I am asking if someone has 
an example of code he wrote to implement a parametric distribution 
from its pdf alone and then used it with MLEstimator().


I am sorry for the post is a bit long but it is a complicate question 
to me and I am not at all skillful in the handling of such notions as 
S4 - class, etc.. so I am a bit lost here..


As a simple example, suppose my theoretical pdf is the skew normal 
distribution (available in package sn):


#skew normal pdf (default values = the standard normal N(0,1)

fsn-function(x,mu=0,sd=1,d=0) {u = (x-mu)/sd;  f = 
dnorm(u)*pnorm(d*u); return(f/sd)}


# d = shape parameter (any real), mu = location (any real), sd = scale 
(positive real)


#to see what it looks like try
x-seq(-1,4,length=200);plot(fsn(x,d=3),type=l)

#Now I tried to create the classes SkewNorm and SkewNormParameter 
copying the example for the binomial

##Class:parameters
setClass(SkewNormParameter,
representation=representation(mu=numeric,sd=numeric,d=numeric),
prototype=prototype(mu=0,sd=1,d=0,name=gettext(Parameter of the Skew 
Normal distribution)),

contains=Parameter
)

##Class: distribution (created using the pdf of the skew normal 
defined above)

setClass(SkewNorm,prototype = prototype(
d = function(x, log = FALSE){fsn(x, mu=0, sd=1,d=0)},
param = new(SkewNormParameter),
.logExact = TRUE,.lowerExact = TRUE),
contains = AbscontDistribution
)

#so far so good but then with
setMethod(mu, SkewNormParameter, function(object) obj...@mu)

#I get the following error message:

 Error in setMethod(mu, SkewNormParameter, function(object) 
obj...@mu) : no existing definition for function mu


I don't understand because to me mu is a parameter not a function... 
maybe that is too complex programming for me and I should switch to 
implementing my likelihood by hand with numerical convolutions and 
optim() etc., but I would like to know how to use distr, so if there 
is anyone who had the same problem and solved it, I would be very 
grateful for the hint !


All the best,
Guillaume





--
Dr. Matthias Kohl
www.stamats.de

__
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] Drawing dendrogram

2009-06-19 Thread Matthias Kohl

maybe ?dendrapply is helpful - in particular see the example section.
There is also some code at http://www.stamats.de/Heatmap.R
(with commentary in German)

hth,
Matthias

Seyit Ali Kayis schrieb:
Dear all, 


I would like to draw a dendrogram and mark some parts/branches (by using segments) 
including their labels. If I draw it without specifying the length of x axix, I am able to do that 
(as in My dendrogram 1 of the following codes). However, if I want to specify the x axix, I am not 
able to draw marking line (by using segments) including labels (as in My dendrogram 2 
of the following codes). Is there a way that I can include my labels into as well? Any help is 
deeply appreciated.

 


Kind Regards

 


Seyit Ali

 

 


Codes:


set.seed(201)
x-matrix(rnorm(100, 0, 1), ncol=20, byrow=TRUE)
myclust-  hclust(dist(cor(x), method='euclidian'), method=ave)
mydend-as.dendrogram(myclust)

par(mfrow=c(1,2), mar=c(5,4,4,3), oma=c(1, 3, 1, 3))

plot(mydend, xlim=c(4, -0.2), horiz = TRUE,  main=My Dendrogram 1)

x0- -0.35
y0- 0.5
x1- 1.8
y1- 0.5
segments(x0, y0, x1, y1, col=red, lty=2)

x0- -0.35
y0- 4.4
x1- 1.8
y1- 4.4
segments(x0, y0, x1, y1, col=red, lty=2)

x0- -0.35
y0- 0.5
x1- -0.35
y1- 4.4
segments(x0, y0, x1, y1, col=red, lty=2)

x0- 1.8
y0- 0.5
x1- 1.8
y1- 4.4
segments(x0, y0, x1, y1, col=red, lty=2)



plot(mydend,  xlim=c(4, 0.5), horiz = TRUE, main=My Dendrogram 2, dLeaf=0.35)#, 
yaxt=n)

x0- 0.3
y0- 0.5
x1- 1.8
y1- 0.5
segments(x0, y0, x1, y1, col=red, lty=2)

x0- 0.3
y0- 4.4
x1- 1.8
y1- 4.4
segments(x0, y0, x1, y1, col=red, lty=2)

x0- 0.3
y0- 0.5
x1- 0.3
y1- 4.4
segments(x0, y0, x1, y1, col=red, lty=2)

x0- 1.8
y0- 0.5
x1- 1.8
y1- 4.4
segments(x0, y0, x1, y1, col=red, lty=2)


-- 
Dr. Seyit Ali KAYIS

Selcuk University
Faculty of Agriculture
Kampus, Konya, TURKEY

s_a_ka...@yahoo.com,s_a_ka...@hotmail.com
Tell: +90 332 223 2830  Mobile: +90 535 587 1139  Fax: +90 332 241 0108

   Greetings from Konya, TURKEY
http://www.ziraat.selcuk.edu.tr/skayis/
-- 








_
[[elided Hotmail spam]]

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


--
Dr. Matthias Kohl
www.stamats.de

__
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] Heatmap: draw horizontal line

2009-05-06 Thread Matthias Kohl

Is

heatmap.2(mat,Rowv=NA,trace = 'none',key=F, add.expr = abline(h = 
c(3.5,5.5), lwd = 3))


what your are looking for?

hth,
Matthias

Paul Evans schrieb:

Hi,

I wanted to draw a heatmap with some horizontal lines. For example:

#-- code --
mat - matrix(-1:1,7,9)
heatmap.2(mat,Rowv=NA,trace = 'none',key=F)

#end code -

In this heatmap, I want to subgroup the rows. For instance, I would like to group rows 5,6 
 7 together, and would draw a black line between rows 4  5 to demarcate the group. 
Similarly, I want to group rows 3  4 in one group, and would like to put a black line 
between rows 2  3.

I guess I'm looking for the 'hline' parameter or something, but couldn't find 
it in the documentation.

Any help would be appreciated!



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


--
Dr. Matthias Kohl
www.stamats.de

__
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] any other fast method for median calculation

2009-04-14 Thread Matthias Kohl
there is function rowMedians in Bioconductor package Biobase which works 
for numeric matrices and might help.

Matthias

Dimitris Rizopoulos wrote:

S Ellison wrote:

Sorting with an appropriate algorithm is nlog(n), so it's very hard to
get the 'exact' median any faster. However, if you can cope with a less
precise median, you could use a binary search between max(x) and min(x)
with low tolerance or comparatively few iterations. In native R, though,
that isn;t going to be fast; interpreter overhead will likely more than
wipe out any reduction in number of comparisons.

In any case, it looks like you are not constrained by the median
algorithm, but by the number of calls. You might do a lot better with
apply, though

apply(df,2,median)


well, for data frames, I think sapply(...) or even unlist(lapply(...)) 
will be faster, e.g.,


mat - matrix(rnorm(50*2e05), 50, 2e05)
DF - as.data.frame(mat)

invisible({gc(); gc()})
system.time(apply(DF, 2, median))

invisible({gc(); gc()})
system.time(sapply(DF, median))

invisible({gc(); gc()})
system.time(unlist(lapply(DF, median), use.names = FALSE))


Best,
Dimitris



On my system 200k columns were processed in negligible time by apply
and I'm still waiting for mapply.

S




Zheng, Xin (NIH) [C] zheng...@mail.nih.gov 14/04/2009 05:29:40


Hi there,

I got a data frame with more than 200k columns. How could I get median
of each column fast? mapply is the fastest function I know for that,
it's not yet satisfied though.
It seems function median in R calculates median by sort and mean.
I am wondering if there is another function with better algorithm.

Any hint?

Thanks,

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


***
This email and any attachments are confidential. Any use...{{dropped:8}}

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





--
Dr. Matthias Kohl
www.stamats.de

__
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 Lambda in Poisson distribution

2009-03-02 Thread Matthias Kohl

Hi,

assuming the assumption of Poisson distribution is true (or at least 
approximately true), you can use MLE as mentioned by Roland using 
(fitdistr of MASS, mle of stats4, MLEstimator of distrMod, probably more).


Using package distrMod you can compute minimum distance estimators (cf. 
function MDEstimator).


Using package ROptEst you can compute (in some sense) optimally robust 
estimators (cf. function roptest).


Best,
Matthias


Rau, Roland wrote:
Hi, 

  

-Original Message-
From: r-help-boun...@r-project.org 
[mailto:r-help-boun...@r-project.org] On Behalf Of Saeed Ahmadi

Sent: Monday, March 02, 2009 3:16 PM
To: r-help@r-project.org
Subject: [R] Finding Lambda in Poisson distribution


Hi,

I have a dataset. First of all, I know that my dataset shall 
follow the

Poission distribution. Now I have two questions:
1) How can I check that my data follow the Poisson distribution?
2) How can I calculate Lambda of my data?



is this maybe some homework?

For 2): I simulated data and did some simple MLE. I don't know if this
is the recommended way, but it worked fine for me. 


Best,
Roland

--
This mail has been sent through the MPI for Demographic Research.  Should you 
receive a mail that is apparently from a MPI user without this text displayed, 
then the address has most likely been faked. If you are uncertain about the 
validity of this message, please check the mail header or ask your system 
administrator for assistance.

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


--
Dr. Matthias Kohl
www.stamats.de

__
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 normalize to a set of internal references

2009-03-01 Thread Matthias Kohl

you should better ask this question on the Bioconductor mailing list.
For qPCR normalisation strategies take a look at
http://www.gene-quantification.info/

Best,
Matthias

Waverley wrote:

Thanks for the advice.  My question is more on how to do this?

Let me use a biology gene analysis example to illustrate:
In biology, there are always some house keeping genes which differ
little even at pathological conditions.

We know that at different batches, there are external factors affect
the measurements.  For example, overall signal intensity might be
different due to lab reagents.
A simplified picture:
Day 1:  Using control samples, I have measured #1 to #110 genes and get data.
Day 2: Using disease samples, I have measured again #1 to #110 genes
and get data.

For those two data sets, I noticed the overall signal intensity in Day
1, for each gene, is more than Day 2.
I know, from biological literature,  gene 101 to 110, are house
keeping genes, should not change much between disease and control.
My questions arise, technically, how do I use gene 101 to 110 values
to adjust the signals of gene 1 to 100 such that the batch effect can
be corrected.  The differences revealing from the comparative analysis
of 1 ~ 100 genes between disease and control are due to biology rather
than lab artifacts.

So the question is how to do that mathematically? If I have only one
house keeping gene, then I can divide every gene to that to normalize,
then compare.  But now I have 10 genes which can be utilized for
normalization.  I assume, the more reference genes to be  used, the
better, under this context.

Can you help again?

Thanks much in advance.


Waverley wrote:
  

Hi,

I have a question of the method as how to normalize the data sets
according to a set of the internal measurements.

For example, I have performed two batches of experiments contrasting
two different conditions (positive versus negative conditions): one at
a time.

1. each experiment, I measure signals of variable v1 to v100. I want
to understand v1 to v100 change under these two contrasting conditions

2. Also I know different variables v101 to v1110, total of 10 of them,
although they are different from each other, but they would of the
same or similar values under these two contrasting conditions

3. How do I do the internal normalization?  How can I use the the
variable v101 to v110 values to normalize the measures of v1 to v100
at either positive or negative condition to minimize batch effect?  I
hope the comparisons of values (v1 to v100) between two different
conditions can be more accurate and robust to external noises.

In general, I have a couple of matrices of the same dimensions and a
reference matrix of values to be used as reference values to be
normalize to.  How should I do that?




I don't understand your problem well, but in general internal
normalization is by and large an attempt to avoid appropriate modeling
(e.g., incorporating block effects or certain covariates in a regression
model), and results in overstated confidence of the final estimates by
not taking into account the imprecision in the normalizing factors.

Frank
  


--
Dr. Matthias Kohl
www.stamats.de

__
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] Does The Randvar package contain a virus(Malware) ?

2009-03-01 Thread Matthias Kohl
there were only minor changes in the latest version of RandVar (from 
0.6.6 to 0.6.7).

Might this be a mirror problem?

Best
Matthias

Timthy Chang wrote:

Today,I update the packages in R.
but AntiVir Guard dectects the  Randvar package as affected file.
What happen ?? 
Thank you for your answer.


http://www.nabble.com/file/p22281513/qq.gif 
  


--
Dr. Matthias Kohl
www.stamats.de

__
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] heatmap.2 color issue

2009-01-21 Thread Matthias Kohl

the heatmapCol function of package MKmisc might help ... try:
library(MKmisc)
example(heatmapCol)

Best
Matthias

Liu, Hao [CNTUS] wrote:

Dear All:

I tried to use heatmap.2 to generate hierarchical clustering using the 
following command:

heatmap.2(datamatrix, scale=row, trace=none, col=greenred(256), 
labRow=genelist[,1], margins=c(10,10), Rowv=TRUE, Colv=TRUE)

datamatrix is subset of a RMA normalized data subset by a genelist. 


The problem is a lot of times, the z-score in key are from, like -5 to 15 or 
-15 to 5, as a result, the zero of z distribution are are either green region 
or red region of
the key, the resulting heatmap are either generally greenish or redish.

I wonder if there is a way to make the heatmap more balanced between red and 
green, I tried to read the heatmap.2 help but could not get a clear idea.

Thanks
Hao


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


--
Dr. Matthias Kohl
www.stamats.de

__
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] Drawing from an empirical distribution

2009-01-06 Thread Matthias Kohl

one could also use package distr; e.g.,

library(distr)
x - 1:10
D - DiscreteDistribution(x)
## = r, d, p and q functions (also with log-argument)
r(D)(5)
p(D)(4)
d(D)(1)
q(D)(0.3)

Best,
Matthias


roger koenker wrote:

Sure, but it would be more 'fun' to modify ecdf() slightly to produce
an ecqf()  function -- essentially reversing the arguments to 
approxfun()--

and then use

ecqf(runif(whatever))

no nit-picking about efficiency, please.


url:www.econ.uiuc.edu/~rogerRoger Koenker
emailrkoen...@uiuc.eduDepartment of Economics
vox: 217-333-4558University of Illinois
fax:   217-244-6678Champaign, IL 61820



On Jan 6, 2009, at 4:42 PM, Antonio, Fabio Di Narzo wrote:


If the ecdf is 'ecdf(x)', do just:

sample(x, size=whatever, replace=TRUE)


HTH,
Antonio.

2009/1/6 culpritNr1 ig2ar-s...@yahoo.co.uk:


Hi All,

Does anybody know if there is a simple way to draw numbers from an 
empirical

distribution?

I know that I can plot the empirical cumulative distribution 
function this

easy:
plot(ecdf(x))

Now I want to pick a number between 0 and 1 and go back to domain of x.
Sounds simple to me.

Any suggestion?

Thank you,

Your culprit
(everybody needs a culprit)


--
View this message in context: 
http://www.nabble.com/Drawing-from-an-empirical-distribution-tp21320810p21320810.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.





--
Antonio, Fabio Di Narzo
Ph.D. student at
Department of Statistical Sciences
University of Bologna, Italy

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


--
Dr. Matthias Kohl
www.stamats.de

__
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] registering a generic method for a class

2009-01-02 Thread Matthias Kohl
In case of S4 objects show is called. Hence, you should implement new 
S4-methods for show ...


Matthias


m.u.r. wrote:

Sorry in advance if this is too simple a question, but I'm stuck with
some odd behavior and I can't find the text to rid myself of this
(admittedly somewhat trivial) problem.  Note that I've done generic
programming with S3 objects in the past, but I've never really
dabbled in creating S4 objects until now.

So, I've created a new S4 class, call it myclass, e.g. :

  

setClass(myclass, representation(slot1 = list));



And I've also created an implementation of print for this class via:

  

setMethod(print, myclass, function(x) { print(paste(a myclass instance with a list of 
length, length(x...@slot1), in slot1)); });



Now I can create a new object, say via:

  

myobj - new(myclass, slot1 = list());



The question now is the difference in output when I type the object
name alone on the command prompt vs when I type print(myobj) on the
command prompt, as seen below:

  

myobj


An object of class “myclass”
Slot slot1:
list()

  

print(myobj)


[1] a myclass instance with a list of length 0 in slot1

I thought that the print method is called when an object is called
from the command line, and to check this I checked the registered
implementations of print via:

  

methods(print)



and in the long output I didn't see the print.myclass version.  I
had thought that the setMethod method would essentially create this
function for me... but I must be missing a key step somewhere along
the line.  Can anyone help me figure out what I would need to do to
get the generic version of the method (in this toy example, print)
to see my class' method as an implementation of the generic?  Thanks
much for any help!

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


--
Dr. Matthias Kohl
www.stamats.de

__
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] Computational Probability

2008-12-26 Thread Matthias Kohl

indeed one could use package distr ...

library(distr)
X - Unif(Min = 0, Max = 1)
Y - convpow(X, 10)
p(Y)(6) - p(Y)(4)

Best
Matthias

Prof Brian Ripley wrote:
Look at packages distr* : they can do your example and might do what 
your real applications.


On Fri, 26 Dec 2008, Rory Winston wrote:


Hi

Firstly , happy Christmas to R-Help! Secondly, I wonder if anyone can 
help

me with the following query: I am trying to reproduce some explicit
probability calculations performed in APPL (a Maple extension for
computational probability). For instance, in APPL, to compute the
probability that the sum of 10 iid uniform variables [0,1] will be 
between 4

and 6, (i..e Pr( 4  \sum_{i=1}^{10}X_i  6)), I can type:

X := UniformRV(0, 1);
Y := ConvolutionIID(X, 10);
CDF(Y,6) - CDF(Y,4);

which gives the required probability .7222. Is there any way to perform
these type of calcuations in R in a general way? I realise that a lot 
of the
machinery behind these computations comes from Maple's symbolic 
engine, but

are there any R extensions for these kind of calculation?

Cheers
Rory

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





--
Dr. Matthias Kohl
www.stamats.de

__
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] Sweave/Rweave and results=verbatim

2008-12-14 Thread Matthias Kohl

it should work without quotes; i.e.,

results = verbatim

hth,
Matthias

Oliver Bandel wrote:

Hello,


in a Rnw-file I have this used stuff to try out tex-results...


==
=
texme - function() cat( {\\bf Hallo, das ist voll fett!}\n )
@

results=verbatim=
texme()
@
==

I used this command: R CMD Sweave example.Rnw


and got this error:


==
[...]
 7 : echo
 8 : echo term verbatim (label=)
Error in match.arg(options$results, c(verbatim, tex, hide)) :
  'arg' should be one of verbatim, tex, hide
Calls: Sweave - SweaveParseOptions - check - match.arg
Execution halted
==


What did I do wrong here?


Any hint is welcome,
   Oliver

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


--
Dr. Matthias Kohl
www.stamats.de

__
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] check if a certain ... argument has been passed on to my user-defined function

2008-12-11 Thread Matthias Kohl

and take a look at
?hasArg

hth,
Matthias

Charles C. Berry wrote:



See

?match.call

and note the expand.dots arg.

HTH,

Chuck


On Thu, 11 Dec 2008, Mark Heckmann wrote:


Hi,

How can I check if a certain ... argument has been passed on to my
user-defined function or not?

foo - function(data, ...)
{
### here I want to check whether xlab was passed with the ... arguments
### or if the ... arguments did not contain an xlab argument
}

I tried  missing(xlab)  , exists(xlab) and several other things but 
did not

find a solution.

TIA,
Mark

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



Charles C. Berry(858) 534-2098
Dept of Family/Preventive 
Medicine

E mailto:cbe...@tajo.ucsd.eduUC San Diego
http://famprevmed.ucsd.edu/faculty/cberry/  La Jolla, San Diego 
92093-0901


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


--
Dr. Matthias Kohl
www.stamats.de

__
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] S4 slot containing either aov or NULL

2008-11-26 Thread Matthias Kohl

Dear Thomas,

take a look at setOldClass ...

## register aov (an S3-class) as a formally defined class
setOldClass(aov)

## list and some others are special cases; cf.
class ? list

## now your code should work
setClassUnion(aovOrNULL, c(aov, NULL))
setClass(c1, representation(value = aovOrNULL))
y1 - new(c1, value = NULL)

#trying to assign an aov object to the slot doesn't work
utils::data(npk, package=MASS)
npk.aov - aov(yield ~ block + N*P*K, npk)
y2 - new(c1, value = npk.aov )

Best,
Matthias

Thomas Kaliwe wrote:

Dear listmembers,

I would like to define a class with a slot that takes either an object 
of class aov or NULL. I have been reading S4 Classes in 15 pages more 
or less and Lecture: S4 classes and methods


#First I tried with list and NULL
setClass(listOrNULL)
setIs(list, listOrNULL)
setIs(NULL, listOrNULL)
#doesn't work

#defining a union class it works with list and null
setClassUnion(listOrNULL, c(list, NULL))
setClass(c1, representation(value = listOrNULL))
y1 = new(c1, value = NULL)
y2 = new(c1, value = list(a = 10))

#but it won't work with aov or null
setClassUnion(aovOrNULL, c(aov, NULL))
setClass(c1, representation(value = aovOrNULL))
y1 = new(c1, value = NULL)

#trying to assign an aov object to the slot doesn't work
utils::data(npk, package=MASS)
npk.aov - aov(yield ~ block + N*P*K, npk)
y2 = new(c1, value = npk.aov )

Any ideas?

Thank you

Thomas Kaliwe

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


--
Dr. Matthias Kohl
www.stamats.de

__
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] Mallows' distance or Earth Mover's distance in R?

2008-10-26 Thread Matthias Kohl

Hi Rainer,

we have not implemented the Mallow's distance so far. Out package 
distrEx includes implementations of the Hellinger, the Kolmogorov, the 
total variation and the Cramer von Mises distance. Combined with our 
package distrMod one can compute minimum distance estimators.


library(distrMod)
?MCEstimator
?MDEstimator

hth,
Matthias

Rainer M Krug wrote:

Hi

I am looking for an implementation (or alternative to) Mallow's
distance or the Earth Mover's distance to compare distributions or
unnormalized distributions (signatures).

Is there an implementation in R or can somebody recommend an alternative?

Thanks

Rainer


  


--
Dr. Matthias Kohl
www.stamats.de

__
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] folded normal distribution in R

2008-10-20 Thread Matthias Kohl

Dear Andreas,

try:
library(distr)
X - Norm() ## standard normal distribution
Y - abs(X)
plot(Y)

Best regards,
Matthias

Andreas Wittmann wrote:

Dear R useRs,

i wanted to ask if the folded normal destribution (Y = abs(X) with X 
normal distributed)
with density and random number generator is implemented in R or in any 
R-related package
so far? Maybe i can use the non-central chi-square distribution and 
rchisq(n, df=1, ncp0) here?


Thanks and best regards

Andreas

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


--
Dr. Matthias Kohl
www.stamats.de

__
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] plot - central limit theorem

2008-10-15 Thread Matthias Kohl

Hi Jörg,

our package distrTeach has functions to visualize the central limit 
theorem and the law of large numbers for arbitrary univariate 
distributions; e.g.,


library(distrTeach)
D - sin(Norm()) + Pois()
plot(D)
illustrateCLT(D, len = 10)
illustrateLLN(D, m = 10)

Best
Matthias

Jörg Groß wrote:

Hi,


Is there a way to simulate a population with R and pull out m samples, 
each with n values

for calculating m means?

I need that kind of data to plot a graphic, demonstrating the central 
limit theorem

and I don't know how to begin.

So, perhaps someone can give me some tips and hints how to start and 
which functions to use.




thanks for any help,
joerg

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


--
Dr. Matthias Kohl
www.stamats.de

__
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 normally distributed values within range

2008-10-06 Thread Matthias Kohl

one could also use the Truncate-methods of package distr; cf.
http://tolstoy.newcastle.edu.au/R/e5/help/08/09/1870.html

If the situation is extreme Peter Daalgard gave some sophisticated code at
http://tolstoy.newcastle.edu.au/R/e5/help/08/09/1892.html

Best,
Matthias

Rolf Turner wrote:


On 7/10/2008, at 11:54 AM, Achaz von Hardenberg wrote:


Hi all,
I need to create 100 normally distributed random values (X) which can
not exceed a specific range (i.e. 0XY).
With rnorm I cannot specify Max and min  values among which values
have to stay, like in runif so does some other simple way exist to do
this with normally distributed random values?


Presumably you want a truncated normal distribution.  Duncan Murdoch
posted some neat code on this list, back in the end of July this year,
to generate samples from such a distribution.


See:

http://finzi.psych.upenn.edu/R/Rhelp02a/archive/137154.html

cheers,

Rolf Turner

##
Attention:\ This e-mail message is privileged and confid...{{dropped:9}}

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


--
Dr. Matthias Kohl
www.stamats.de

__
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] help on sampling from the truncated normal/gamma distribution on the far end (probability is very low)

2008-09-18 Thread Matthias Kohl

you could use package distr and function Truncate; e.g.

library(distr)
N - Norm(mean = -4, sd = 1)
NT - Truncate(N, lower = 0, upper = Inf)
r(NT)(10)

Unfortunatelly, your example using sd = 0.1 is very extreme and Truncate 
doesn't work; see also
pnorm(0, mean = -4, sd = 0.1, lower.tail = FALSE) == 0 ## which on my 
system is TRUE


Best,
Matthias

Moshe Olshansky wrote:

Well, I made a mistake - your lambda should be 400 and not 40!!!


--- On Thu, 18/9/08, Moshe Olshansky [EMAIL PROTECTED] wrote:

  

From: Moshe Olshansky [EMAIL PROTECTED]
Subject: Re: [R] help on sampling from the truncated normal/gamma distribution 
on the far end (probability is very low)
To: r-help@r-project.org, Daniel Davis [EMAIL PROTECTED]
Received: Thursday, 18 September, 2008, 5:00 PM
Hi Sonia,

If I did not make a mistake, the conditional distribution
of X given that X  0 is very close to exponential
distribution with parameter lambda = 40, so you can sample
from this distribution.


--- On Mon, 15/9/08, Daniel Davis
[EMAIL PROTECTED] wrote:



From: Daniel Davis [EMAIL PROTECTED]
Subject: [R] help on sampling from the truncated
  

normal/gamma distribution on the far end (probability is
very low)


To: r-help@r-project.org
Received: Monday, 15 September, 2008, 2:28 PM
Hi, guys,

I am trying to sample from a truncated normal/gamma
distribution.
But only the far end of the distribution (where the
probability is very low)
is left. e.g.

mu = - 4;
sigma = 0.1;
The distribution is Normal(mu,sigma^2) truncated on
[0,+Inf];

How can I get a sample? I tried to use inverse CDF
  

method,


but got Inf as
answers. Please help me out.

Also, pls help me on the similar situation on gamma
dist'n.


Thanks,
Sonia

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


--
Dr. Matthias Kohl
www.stamats.de

__
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] All possible pairs of two variables

2008-09-11 Thread Matthias Kohl
take a look at the code of the function pairwise.table. It shows a way 
using function outer.


hth,
Matthias


Yihui Xie wrote:

Just simple loops as follows:


z = NULL
for (y in 0:10) {
for (x in 0:y) {
z = rbind(z, c(x, y))
}
}
z


Yihui

On Thu, Sep 11, 2008 at 12:34 PM, Ron Michael [EMAIL PROTECTED] wrote:
  

I have two variables (x,y) :

x : it takes all integer values from 0 to y and,
y : takes all values from 0, 10

I am looking for some R code to find all possible pairs of (x,y). Can anyone 
please help me?






  


--
Dr. Matthias Kohl
www.stamats.de

__
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] fft: from characteristic funtion to density function

2008-09-11 Thread Matthias Kohl

take a look at
http://finzi.psych.upenn.edu/R/Rhelp02a/archive/130153.html

hth,
Matthias


Jindan Zhou wrote:

Hello all!
I've posted the question before but I am still struggling for an
answer, please help if you can;-)

Suppose a discrete series of data is generated by the following
equation:
CF=exp(-(t^2)/2)
which is the characteristic function of a random variable X with
standard normal distribution, how do I *numerically* solve for its
probability density function? i.e., obtain a series of data that plots
a well know bell-shape curve? Or, find x such that Pr(Xx)=0.05?

I have no problem to derive the forward and inverse Fourier equation
pair, but just no clue on numerical method. A working code on this
simplest exercise
would clear off many, many of my confusions
I have generated CF with t of 256 points equally spaced between -4 to
4 in Excel, but don't know
how to proceed.

Thanks for your help!
Jindan

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


--
Dr. Matthias Kohl
www.stamats.de

__
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] Package for Tidying-up R-code

2008-09-05 Thread Matthias Kohl

http://cran.r-project.org/doc/manuals/R-exts.html#Tidying-and-profiling-R-code

hth,
Matthias

Don MacQueen wrote:
I don't know what kind of tidying you want, but ESS within emacs is a 
possibility.


-Don

At 6:32 PM +0900 9/5/08, Gundala Viswanath wrote:

Dear all,

Is there any such package?

I am thinking of something equivalent to Perl::Tidy.
For example with VI editor one can visual highlight a portion
of Perlcode and then issue the command:

!perltidy

then the code will be automatically arranged.

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





--
Dr. Matthias Kohl
www.stamats.de

__
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] Plotting horizontal dendrograms

2008-06-02 Thread Matthias Kohl

Dear Charles,

the sample code you cite has an additional as.dendrogram compared to 
your code. Hence, dend1 is of class dendrogram (class(dend1)) and if you 
call plot(dend1) the method for dendrograms (i.e., plot.dendrogram) is 
executed which has a horiz argument; confer ?plot.dendrogram.
Whereas in your code h is of class hclust (class(hclust)) and therefore 
plot(h) leads to a call of plot.hclust which has no horiz argument; cf. 
?plot.hclust.
So, adding as.dendrogram in your code should remove the warning messages 
and give what you want; e.g.


plot(as.dendrogram(h), horiz = TRUE)

Best,
Matthias

Charles K. Minns wrote:

I am using hclust and plot to produce dendrograms. Using my input data I am
able to complete an analysis and obtain a vertical plot.
I want to be able to plot the dendrogram horizontally.I am using version 2.6
of R and have updated my packages recently.
 
Using the sample script for dendrograms I can produce a horizontal plot

using the instruction horiz = TRUE in plot().
When I use the same instruction in my own script I get a warning message
that horiz, and horizontal, are not recognized commands in plot.
I cannot find any documentation in R on how to change the orientation of
plots and general searches for information on this subject.
Can someone explain why my script does not work.
 
Sample code follows:

---
# Sample code extracted from the Dendrogram documentation in R
#
require(graphics); require(utils)
 
hc - hclust(dist(USArrests), ave)

(dend1 - as.dendrogram(hc)) # print() method
str(dend1)  # str() method
str(dend1, max = 2) # only the first two sub-levels
 
op - par(mfrow= c(2,2), mar = c(5,2,1,4))

plot(dend1, horiz = TRUE)
plot(dend1)
#
# My code
#
x - read.table(Rclust2.data,header=TRUE)
 
m - data.matrix(x)

d - as.dist(m)
h - hclust(d)
plot(h, horiz = TRUE)

--
The three plot commands execute. Plot 1 is horizontal, plot 2 is vertical
and plot 3 is vertical. Plus plot 3
generates the console output:
-
  

#
# My code
#
x - read.table(Rclust2.data,header=TRUE)

m - data.matrix(x)
d - as.dist(m)
h - hclust(d)
plot(h, horiz = TRUE)


Warning messages:
1: In plot.hclust(h, horiz = TRUE) : horiz is not a graphical parameter
2: In plot.hclust(h, horiz = TRUE) : horiz is not a graphical parameter
3: In title(main = main, sub = sub, xlab = xlab, ylab = ylab, ...) :
  horiz is not a graphical parameter
---
 
I would appreciate help solving this problem. Doubtless there is a simple

answer.

 


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


--
Dr. Matthias Kohl
www.stamats.de

__
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] Eliminating [...] from print

2008-06-01 Thread Matthias Kohl

Dear Ted,

what about
cat(\n)

or

message()

Best,
Matthias

(Ted Harding) wrote:

Hi Folks,
This must be easy but I've not managed to locate the solution!

Basically: I'm using sink() to save successively obtained
results, e.g. I construct a set of regression coefficients
etc. with names rows and columne (I want to see the names
in the output) as an object (say Object), and then I
emit it with

  print(Object)

So far so good. But I also want to print a little header
and separator prior to each such print(Object) which would
identify it and explain what's going on. BUT: I can't get
rid of the [1] which comes out each time I print a
character string. Getting rid of the quotes is easy: just
use quote=FALSE. For example:

  print (,quote=FALSE)

gives

  [1] 

But how do I get rid of that [1]??

With thanks,
Ted.



E-Mail: (Ted Harding) [EMAIL PROTECTED]
Fax-to-email: +44 (0)870 094 0861
Date: 01-Jun-08   Time: 18:53:00
-- XFMail --

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


--
Dr. Matthias Kohl
www.stamats.de

__
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] Heatmap.2 - eliminate cluster and dendrogram

2008-05-15 Thread Matthias Kohl

Dear Peter,

heatmap.2(z, dendrogram = none, Rowv = FALSE, Colv = FALSE)

should do what you want.

Best,
Matthias

Peter Scacheri wrote:
Using the heatmap.2 function, I am trying to generate a heatmap of a 2 
column x 500 row matrix of numeric values.  I would like the 1st 
column of the matrix sorted from the highest to the lowest values - so 
that the colors reflected in the first column of the heatmap (top to 
bottom) go from red to green.


After sorting the matrix (z), I tried the following command, but the 
data remains clustered.


heatmap.2((z),col=greenred(100),dendrogram=NULL,Rowv=FALSE)

I also tried the following, but the data remained clustered.

heatmap.2((z),col=greenred(100),Colv=1:ncol(z))

ANY IDEAS??
Thanks so much,
Peter

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


--
Dr. Matthias Kohl
www.stamats.de

__
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] Function for subset of cases/lines

2008-05-15 Thread Matthias Kohl

Dear Stefan,

what about

q1-c(4660,5621,5629,8030,8080,8180,8501,8190,8370,8200)
x - data.frame(number = c(4660, 5010, 5621, 5629, 8030, 8080 , 8090, 
8180, 8501, 8190, 8200, 8370, 8200), height = c(2.5, 1.4, 0.8, 2.3, 2.5, 
2.4, 0.9, 1.4, 1.2, 1.9, 2.0, 2.1, 1.8))


ind - x[,number] %in% q1
mean(x[ind,height])

## or shorter
mean(x[x[,number] %in% q1,height])

Best,
Matthias

Stefan Uhmann wrote:

Hi,

I have a vector:

q1-c(4660,5621,5629,8030,8080,8180,8501,8190,8370,8200)

The following command gives me the mean of its elements:

mean(q1)
[1] 7346.1

What can I do to do the same for the variable 'height', but only for 
the cases/rows which have one of the elements of q1 as 'number':


  number height
1 46602.5
2 50101.4
3 56210.8
4 56292.3
5 80302.5
6 80802.4
7 80900.9
8 81801.4
9 85011.2
10 81901.9
11 82002.0
12 83702.1
13 82001.8

Thanks for any help,
Stefan

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


--
Dr. Matthias Kohl
www.stamats.de

__
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] Heatmap.2 - eliminate cluster and dendrogram

2008-05-15 Thread Matthias Kohl

Dear Peter,

the call works for me; e.g.

x - matrix(rnorm(100), ncol = 4)
heatmap.2(x, dendrogram=none, Colv = FALSE, Rowv = FALSE)

There seems to be something wrong with your matrix z ...

Best,
Matthias

Peter Scacheri wrote:

Thanks Matthias,

I tried that, but received the following error message:

Error in image.default(1:nc, 1:nr, x, xlim = 0.5 + c(0, nc), ylim = 
0.5 +  :

dimensions of z are not length(x)(+1) times length(y)(+1)


-Peter

At 10:18 AM +0200 5/15/08, Matthias Kohl wrote:

Dear Peter,

heatmap.2(z, dendrogram = none, Rowv = FALSE, Colv = FALSE)

should do what you want.

Best,
Matthias

Peter Scacheri wrote:
Using the heatmap.2 function, I am trying to generate a heatmap of a 
2 column x 500 row matrix of numeric values.  I would like the 1st 
column of the matrix sorted from the highest to the lowest values - 
so that the colors reflected in the first column of the heatmap (top 
to bottom) go from red to green.


After sorting the matrix (z), I tried the following command, but the 
data remains clustered.


heatmap.2((z),col=greenred(100),dendrogram=NULL,Rowv=FALSE)

I also tried the following, but the data remained clustered.

heatmap.2((z),col=greenred(100),Colv=1:ncol(z))

ANY IDEAS??
Thanks so much,
Peter

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


--
Dr. Matthias Kohl
www.stamats.de


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


--
Dr. Matthias Kohl
www.stamats.de

__
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] dataframes to a list

2008-05-14 Thread Matthias Kohl

Does

mget(paste(data, 1:50, sep = ), envir = globalenv())

give what you want?

hth,
Matthias

Shubha Vishwanath Karanth wrote:

How do I proceed from this stage?
paste(data,1:10,sep=,collapse=,)

Shubha Karanth | Amba Research
Ph +91 80 3980 8031 | Mob +91 94 4886 4510 
Bangalore * Colombo * London * New York * San José * Singapore * www.ambaresearch.com


-Original Message-
From: [EMAIL PROTECTED] [mailto:[EMAIL PROTECTED] On Behalf Of Shubha 
Vishwanath Karanth
Sent: Wednesday, May 14, 2008 3:59 PM
To: [EMAIL PROTECTED]
Subject: [R] dataframes to a list

 


Hi R,

 


I have the data frames, data1, data2data50. Now I want to put all of
these in a single list. But,

 

  

list(data1, data2,.data50) is very big to write. How do I then


do it?

 


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.
This e-mail may contain confidential and/or privileged i...{{dropped:10}}

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


--
Dr. Matthias Kohl
www.stamats.de

__
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] call function immediately before plot.new()

2008-04-28 Thread Matthias Kohl

Hi Jake,

are you looking for argument panel.first of plot.default?
?plot.default

panel.first: an expression to be evaluated after the plot axes are set
 up but before any plotting takes place.  This can be useful
 for drawing background grids or scatterplot smooths.

hth,
Matthias


Jake Michaelson wrote:

Hi all,

I would like to be able to call a custom function automatically before
plot.new() is called (more specifically, before a new plot is created on the
current graphics device).  Recently I've been poking around in the help
files of some of the low(er) level plotting functions, and I seem to
remember something saying that you could somehow add a setting to call a
function immediately before plot.new() is called (something kind of in the
spirit of .First).  I've been looking everywhere, and I can't find the
documentation I'm thinking of, which makes me wonder if I'm just making all
this up.

Can this be done?

Thanks,

Jake

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


--
Dr. Matthias Kohl
www.stamats.de

__
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] Sum of random values

2008-04-27 Thread Matthias Kohl
Hello Ricardo,

no, the distribution needs not to be predefined. You can define your own 
distribution. For more details you may want to have a look at the 
vignette which is provided by package distrDoc.

Do the random variables follow a discrete or an absolutely continuous 
distribution?
To see how convolution is performed in package distr have a look at

library(distr)
getMethod(+, c(DiscreteDistribution, DiscreteDistribution))
## respectively
getMethod(+, c(AbscontDistribution, AbscontDistribution))

hth,
Matthias

Ricardo Bessa wrote:
 Thank you Gregory,
 I know the package distR but the the percentiles that I have are from an 
 unkonwn distribuition, and I think that the package distR only works with 
 pre-defined distribuions.
  
 Best Regards,
 Ricardo Bessa Subject: RE: [R] Sum of random values Date: Fri, 25 Apr 2008 
 09:16:10 -0600 From: [EMAIL PROTECTED] To: [EMAIL PROTECTED]; 
 r-help@r-project.org  You might want to look at the distr package (and its 
 relatives). It provides methods for calculating the distribution function 
 of combinations (the sum is one) of other distributions. I'm not sure how 
 you would convert your percentiles to a distribution function, but there may 
 be a way in the documentation for distr and friends.  --  Gregory (Greg) 
 L. Snow Ph.D. Statistical Data Center Intermountain Healthcare [EMAIL 
 PROTECTED] (801) 408-8111 -Original Message-  From: 
 [EMAIL PROTECTED]   [mailto:[EMAIL PROTECTED] On Behalf Of Ricardo Bessa  
 Sent: Thursday, April 24, 2008 5:12 PM  To: r-help@r-project.org  
 Subject: [R] Sum of random values  Hello,  I have two random 
 variables with their percentiles which   correspond to their !
 pr!
  obability distribution function. My   objective is to sum these two random 
 variables. There exists   any algorithm or procedure in R capable of 
 converting the   percentiles to a probability density function? is the fast 
   Fourier transform function of R(fft) capable of doing the sum   with a 
 convolution?I'm just starting with this specific problem, so any help 
 it   will be very useful.Best regards,  Ricardo Bessa
 _
 [[elided Hotmail spam]][[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.   
 _

 [[elided Hotmail spam]]

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

-- 
Dr. Matthias Kohl
www.stamats.de

__
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] inverse F distribution in R?

2008-04-24 Thread Matthias Kohl
just a comment to fstat...

Using package distrEx one also has fstat functionality

library(distrEx)
E(Fd(df1 = 3, df2 = 3))
E(Fd(df1 = 4, df2 = 4))
E(Fd(df1 = 5, df2 = 5))
var(Fd(df1 = 5, df2 = 5))

for more functionals (skewness, kurtosis, IQR, ...) see for instance
help(var, package = distrEx)

Best
Matthias

Peter Dalgaard wrote:
 Jennifer Balch wrote:
   
 Dear all,

 I'm looking for a function that calls the inverse F-distribution.  
 Something equivalent to FINV in matlab or excel.

 Does anyone know if such a function already exists for R? (I haven't  
 been able to find one.)

 Thanks for any leads.

   
 
 I would guess that they are qf() or a variant thereof. I gather that
 _fpdf
 http://www.mathworks.com/access/helpdesk/help/toolbox/stats/fpdf.html__,
 fcdf
 http://www.mathworks.com/access/helpdesk/help/toolbox/stats/fcdf.html,
 finv
 http://www.mathworks.com/access/helpdesk/help/toolbox/stats/finv.html,
 frnd
 http://www.mathworks.com/access/helpdesk/help/toolbox/stats/frnd.html
 i_n matlab is df, pf, qf, and rf in R. We don't have an equivalent of
 fstat for mean and variance of the F distribution (this might actually
 be nice to have).




   
 Best,
 Jennifer

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


   

-- 
Dr. Matthias Kohl
www.stamats.de

__
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] mean for each row

2008-03-17 Thread Matthias Kohl
or rowMeans ...
Matthias

Hesen Peng wrote:
 How about rowSums(x)/ncol(x), where x is the matrix?


 On Mon, Mar 17, 2008 at 1:48 PM, Roslina Zakaria [EMAIL PROTECTED] wrote:
   
   Hi r-users,
  How do find the mean for each row? Thank you in advance for your help.


   1981 1982 1983 1984 1985 1986 1987 1988 1989 1990 Day Totals
  10.0  0.0  0.0  0.0  0.6  0.0  8.4  0.0 29.4  0.0   38.4
  20.0  0.0  1.8  0.0 22.4  0.0  0.2  0.4  0.8  0.0   25.6
  37.8  0.0  0.0 17.6  1.4  0.0  0.0  0.0  0.0  0.0   26.8
  42.2  0.8  0.4  0.0  0.2 11.2  1.4 33.2  0.0  0.0   49.4
  50.2  1.8  0.0  1.0  0.0  0.2  0.0 12.2  0.0 19.2   34.6
  60.0  0.0  0.0  1.0  0.0  0.0  0.0  2.2  0.0 14.6   17.8
  70.0  0.0  0.0  0.0  3.6  0.2  0.0  2.0  0.0  0.26.0
  80.0  0.0 10.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0   10.0
  90.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.00.0
  10   1.0  1.4  0.0  0.0  0.0  1.0  0.4  0.0  0.0  0.03.8
  11   0.0  8.2  0.0  0.0  0.0  0.0  8.4  0.0  0.0  0.4   17.0
  12  10.8  0.8  0.0  0.0  0.0  0.0  1.0  0.0  0.0  4.2   16.8
  13  32.8  0.0  0.0  0.8  0.0  0.0  0.0  0.0  0.2  0.2   34.0
  14   1.0  0.0  1.6  0.2  0.0  0.0  0.0  0.0  0.0  0.02.8
  15   0.0  0.0  0.0  2.2  1.4  0.0  0.0  0.0  0.0  0.03.6
  16   0.6  0.0  0.0  0.6 22.0  0.0  0.0  0.0  0.0  0.0   23.2
  17   0.0  0.0  0.0  0.0  0.0  0.0  2.8  0.0  0.0  0.02.8
  18   2.8  0.0  0.0  0.0  0.0  0.0  8.2  0.0  0.0  8.2   19.2
  19   0.0  0.0  0.0  0.0  0.0  0.2  0.0  0.0  0.0  0.20.4
  20   0.0  8.2  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.08.2
  21   0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.00.0
  22   0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.40.4
  23   0.0  0.0  0.0  0.0  0.0  0.0  0.0  1.2  0.0  0.01.2
  24   0.0  0.0  0.0  0.0  1.0  0.0  0.0  2.0  8.2  0.0   11.2
  25   3.2  0.0  0.0  0.6  0.6  0.0  0.0  0.0 11.8  0.0   16.2
  26   0.0  0.0 26.2  0.0 12.6  0.0  0.0  2.2  0.0  0.0   41.0
  27   0.2  0.0 10.6  0.0  1.2  0.0  0.0  1.8  0.0  0.0   13.8
  28   0.0  4.0  0.0  5.8  0.0  0.0  0.0  0.0  0.0  0.09.8
  29   0.2 12.4  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0   12.6
  30   0.0  2.6  0.0  0.0  2.2  0.0  0.0  0.0  0.0  0.04.8
  31   0.0  0.0  0.6  0.0  0.8  0.0  0.0  0.0  4.8  0.06.2
  Year Totals 62.8 40.2 51.2 29.8 70.0 12.8 30.8 57.2 55.2 47.6  457.6


   
 
  Be a better friend, newshound, and

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

 



   

-- 
Dr. Matthias Kohl
www.stamats.de

__
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 boxplots using boxplot

2008-03-01 Thread Matthias Kohl
Hi Marek,

use a list instead of a data.frame; e.g.

boxplot(list(x = rnorm(10), y = rnorm(20)))

hth,
Matthias

Marek Bartkuhn wrote:
 Hi,

 I have data in 3 vectors a,b and c looking like

 eg  a:1.2 3.4 1.4 ..


 I like to have the data from each vector as a boxplot, but the 3  
 boxplots within one graph.

 One example is the first graph  at

 http://www.statmethods.net/graphs/boxplot.html

 Unfortunately I do not really understand the corresponding  
 description. The vectors are of different length, therefore R  
 complains when I try to build a data.frame of the 3 vectors:

 arguments imply differing number of rows: 6175, 20, 190

 What do I have to do?

 Marek

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

-- 
Dr. Matthias Kohl
www.stamats.de

__
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] likelihood from test result

2008-01-10 Thread Matthias Kohl
Dear David,

no, distrTEst won't help. It has a different intention.

We are currently working on a new package distrMod
(cf. https://r-forge.r-project.org/projects/distrmod/)
which sometime might have such a functionality.

Best,
Matthias

David Bickel wrote:
 Is there any automatic mechanism for extracting a likelihood or test
 statistic distribution (PDF or CDF) from an object of class htest or
 from another object of a general class encoding a hypothesis test
 result?

 I would like to have a function that takes x, an object of class
 htest, as its only argument and that returns the likelihood or test
 statistic distribution that was used to compute the p-value. It seems
 the only way to write such a function is to manually assign each test
 its statistic's distribution, e.g., like this:

 FUN - if(names(x$statistic) == t)
   dt
 else if(names(x$statistic) == X-squared)
   dchisq
 # etc.

 Is there a general S3 or S4 class other than htest that would better
 accommodate such extraction of distributions or likelihoods? I would
 also appreciate any suggestions for strategies or contributed packages
 that may facilitate automation. For example, would the distrTEst
 package help?

 David

 __
 David R. Bickel
 Ottawa Institute of Systems Biology
 http://www.oisb.ca/members.htm

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


__
R-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] probability from different values

2007-12-27 Thread Matthias Kohl
one addition ...
 Hello Ricardo,

 another solution could be using package distr:

 library(distr)
 A - c(18,18,18,19,20,21,22,23,24,25,26,27,28)
 DA - DiscreteDistribution(A)
 # maybe
 # support(DA)
 # plot (DA)

 B - c(82,83,84,85,85,86,87,88,89,90,91,91,92)
 DB - DiscreteDistribution(B)
 # support(DB)
 # plot(DB)

 DC - DB - DA# convolution with of DB with (-DA)
 # support(DC)
 # plot(DC)
   
to get the probabilities for the C-values:
p(DC)(support(DC))

 hth,
 Matthias

 [EMAIL PROTECTED] wrote:
   
 i have two vectors (A and B) and i need create another vector (C) 
 from the subtraction of A's values with the B's values. How can i 
 estimate the probability of C's values if i have differents values 
 combinations of A and B that can result in the same value? like 90 
 -20 = 70 and 91 - 21 = 70.

 example:
 B = {18  1818 19  20  21  22  23  24
 25  26  27   28} 
 A = {82 83 84 85  85 86   87 88   89
 90  91 9192 }
 
   
 Try this:
 A = c(18,18,18,19,20,21,22,23,24,25,26,27,28)
 B = c(82,83,84,85,85,86,87,88,89,90,91,91,92)
 lenA = length(A)
 lenB = length(B)
 AA = rep(A, lenB)
 BB = rep(B, each=lenA)
 table(AA-BB)/(lenA*lenB)

 Regards,
 Richie.

 Mathematical Sciences Unit
 HSL


 
 ATTENTION:

 This message contains privileged and confidential inform...{{dropped:20}}

 __
 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] probability from different values

2007-12-27 Thread Matthias Kohl
Hello Ricardo,

another solution could be using package distr:

library(distr)
A - c(18,18,18,19,20,21,22,23,24,25,26,27,28)
DA - DiscreteDistribution(A)
# maybe
# support(DA)
# plot (DA)

B - c(82,83,84,85,85,86,87,88,89,90,91,91,92)
DB - DiscreteDistribution(B)
# support(DB)
# plot(DB)

DC - DB - DA# convolution with of DB with (-DA)
# support(DC)
# plot(DC)

hth,
Matthias

[EMAIL PROTECTED] wrote:
 i have two vectors (A and B) and i need create another vector (C) 
 from the subtraction of A's values with the B's values. How can i 
 estimate the probability of C's values if i have differents values 
 combinations of A and B that can result in the same value? like 90 
 -20 = 70 and 91 - 21 = 70.

 example:
 B = {18  1818 19  20  21  22  23  24
 25  26  27   28} 
 A = {82 83 84 85  85 86   87 88   89
 90  91 9192 }
 

 Try this:
 A = c(18,18,18,19,20,21,22,23,24,25,26,27,28)
 B = c(82,83,84,85,85,86,87,88,89,90,91,91,92)
 lenA = length(A)
 lenB = length(B)
 AA = rep(A, lenB)
 BB = rep(B, each=lenA)
 table(AA-BB)/(lenA*lenB)

 Regards,
 Richie.

 Mathematical Sciences Unit
 HSL


 
 ATTENTION:

 This message contains privileged and confidential inform...{{dropped:20}}

 __
 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] Simulate an AR(1) process via distributions? (without specifying a model specification)

2007-11-28 Thread Matthias Kohl
Dear Pedro,

you might be interested in the demo StationaryRegressorDistr of 
package distr.

library(distr)
demo(StationaryRegressorDistr)

hth,
Matthias


[EMAIL PROTECTED] wrote:
 Thanks Prof. Ripley.

 My apologies for not including the code.

 Below I illustrate my point using the GLD package. 

 Thank you very much for your time.

 Kind Regards,

 Pedro N. Rodriguez 


 # Code begins

 # Simulate an ar(1) process  
 # x = 0.05 + 0.64*x(t-1) + e

 # Create the vector x
   x   - vector(length=1000)

 #simulate the own risk
   e   - rnorm(1000)

 #Set the coefficient
   beta- 1.50

 # set an initial value
   x[1]- 5

 #Fill the vector x
   for(i in 2:length(x))
   {
   x[i]- 0.05 + beta*x[i-1] + e[i] 
   }

 #Check the AR(1) 
   simulated_data_ar - arima(x,order=c(1,0,0))
   simulated_data_ar

 #Using the G Lambda Distribution to fit the distribution. 
   library(gld)
   resul1  - starship(x,optim.method=Nelder-Mead) 
   lambdas1- resul1$lambda

 #Plot the Distribution
   plotgld(lambdas1[1],lambdas1[2],lambdas1[3],lambdas1[4])

 #Random Deviates from GLD
   x_sim   -
 rgl(1000,lambdas1[1],lambdas1[2],lambdas1[3],lambdas1[4])

 #Fit an AR(1)
   gld_simulated   - arima(x_sim,order=c(1,0,0))
   gld_simulated

 #Code ends


 -Original Message-
 From: Prof Brian Ripley [mailto:[EMAIL PROTECTED] 
 Sent: Wednesday, November 28, 2007 11:37 AM
 To: Rodriguez, Pedro
 Cc: [EMAIL PROTECTED]
 Subject: Re: [R] Simulate an AR(1) process via distributions? (without
 specifying a model specification)

 On Wed, 28 Nov 2007, [EMAIL PROTECTED] wrote:

   
 Is it possible to simulate an AR(1) process via a distribution?
 

 Any distribution *of errors*, yes.  Of the process values, not in
 general.

   
 I have simulated an AR(1) process the usual way (that is, using a
 
 model
   
 specification and using the random deviates in the error), and used
 
 the
   
 generated time series to estimate 3- and 4-parameter distributions
 
 (for
   
 instance, GLD). However, the random deviates generated from these
 distributions do not follow the specified AR process.
 

 How do you know that?  Please give us the reproducible example we asked 
 for (in the posting guide, at the bottom of every message), and we
 should 
 be able to explain it to you.



__
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] skewness: generating distribution and estimating M3.

2007-11-20 Thread Matthias Kohl
Dear Milton,

package distrEx has methods to compute skewness.

library(distrEx)
?skewness

To generate new distributions you could use package distr.

hth,
Matthias


Milton Cezar Ribeiro wrote:
 Dear all,

 How can I generate values which distrituions have left or rigth skewness.
 After that, which function compute the skewness of each distributions?

 Kind regards.

 miltinho



  para armazenamento!

   [[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] analytical solution to Sum of binominal distributed random numbers?

2007-10-24 Thread Matthias Kohl
Hello Rainer,

You could also get a very exact approximation using our package distr; e.g.

library(distr)
B1 - Binom(prob = 0.4, size = 20)
B2 - Binom(prob = 0.2, size = 10)
B3 - B1+B2 #convolution!
plot(B3)

hth,
Matthias

Rainer M. Krug schrieb:
 Frede Aakmann Tøgersen wrote:
   
 Perhaps

 http://stinet.dtic.mil/cgi-bin/GetTRDoc?AD=ADA266969Location=U2doc=GetTRDoc.pdf

 is something that you can use?
 

 Thanks a lot - that might help.

 Rainer

   

 Best regards

 Frede Aakmann Tøgersen
 Scientist


 UNIVERSITY OF AARHUS
 Faculty of Agricultural Sciences
 Dept. of Genetics and Biotechnology
 Blichers Allé 20, P.O. BOX 50
 DK-8830 Tjele

 Phone:   +45 8999 1900
 Direct:  +45 8999 1878

 E-mail:  [EMAIL PROTECTED]
 Web:http://www.agrsci.org

 This email may contain information that is confidential.
 Any use or publication of this email without written permission from Faculty 
 of Agricultural Sciences is not allowed.
 If you are not the intended recipient, please notify Faculty of Agricultural 
 Sciences immediately and delete this email.

  

  

 
 -Oprindelig meddelelse-
 Fra: [EMAIL PROTECTED] 
 [mailto:[EMAIL PROTECTED] På vegne af Rainer M Krug
 Sendt: 24. oktober 2007 09:11
 Til: Charles C. Berry
 Cc: r-help
 Emne: Re: [R] analytical solution to Sum of binominal 
 distributed random numbers?

 Hi Charles

 thanks for the pointing out that size and prob can be vectors 
 as well - I tried it out but used 1 as the number of 
 observations, assuming that and it only gave me one randon 
 mumbewr (as it should be but not expected).

 But I was more looking at a analytical solution, as I have to 
 sum up a huge number of random numbers. But I am going to try 
 your solution as it should be much faster already.

 Thanks

 Rainer


 Charles C. Berry wrote:
   
?rbinom

 only says:

  size: number of trials (zero or more).

  prob: probability of success on each trial.


 But they can be vectors.

 BTW, you were aked to PLEASE ... provide  minimal, self-contained, 
 reproducible code.

 What you show cannot run without correction.

 Most likely, you intended size(n) to be the n-th element of 
 
 the vector 
   
 'size', which in R is written 'size[ n ]' .

 In which case

sum (rbinom( length(prob) , size, prob ) )

 works.

 Chuck

 On Tue, 23 Oct 2007, Rainer M Krug wrote:

 
 Hi

 I have two vectors, prob and size, and I want to add the random 
 deviates  of these two, i.e.

 sum(
   sapply(
  1:length(prob),
  function(n){ rbinom(1, size(n), prob(n) }
 )
 )

 My problem is that I have to do this for a large number of value 
 combinations. Is there a faster way of doing this?

 Rainer

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

   
 Charles C. Berry(858) 534-2098
  Dept of 
 
 Family/Preventive Medicine
   
 E mailto:[EMAIL PROTECTED] UC San Diego
 http://famprevmed.ucsd.edu/faculty/cberry/  La Jolla, San Diego 
 92093-0901

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