Re: [R] Variable name as string

2010-10-17 Thread Jan private
So here is the next version.

Why does the intercept needs lower.tail=TRUE to give the same result as
summary() for value=0?

# See Verzani, simpleR (pdf), p. 80
coeff.test - function(lm.result, idx, value) {
  # idx = 1 is the intercept, idx1 the other coefficients
  # null hypothesis: coeff = value
  # alternative hypothesis: coeff != value
  coeff - coefficients(lm.result)[idx]
  SE - coefficients(summary(lm.result))[idx,Std. Error]
  n - df.residual(lm.result)
  t - (coeff - value )/SE
  if (idx == 1) {
2 * pt(t,n,lower.tail=TRUE) # times two because problem is two-sided
  } else {
2 * pt(t,n,lower.tail=FALSE)
  }
}

__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/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] Variable name as string

2010-10-16 Thread Jan private
Hello,

from Verzani, simpleR (pdf), p. 80, I created the following function to
test the coefficient of lm() against an arbitrary value.

coeff.test - function(lm.result, var, coeffname, value) {
  # null hypothesis: coeff = value
  # alternative hypothesis: coeff != value
  es - resid(lm.result)
  coeff - (coefficients(lm.result))[[coeffname]]
  # degrees of freedom = length(var) - number of coefficients?
  n - df.residual(lm.result) 
  s - sqrt( sum( es^2 ) / n )
  SE - s/sqrt(sum((var - mean(var))^2))
  t - (coeff - value )/SE
  2 * pt(t,n,lower.tail=FALSE) # times two because problem is two-sided
}

E.g. it needs to be called coeff.test(lm(N ~ D), D, D, 70) if I want
to test the probability slope of the linear model being 70.

Is there a more elegant way to avoid passing D and D both as
parameters?

Also, as a non-professional, I would like to know whether the function
is valid for all coefficients of lm(), e.g. coeff.test(lm(N ~ D + H), H,
H, 70). I am aware that Verzani gives a different formula for testing
the intercept.

Thanks!
Jan Rheinländer

__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/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] Variable name as string

2010-10-16 Thread Jan private

 coeffname - deparse(substitute(var))
 
 Is that what you wanted?
 
Yes! That gets rid of the extra parameter.

Thanks

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


Re: [R] Odp: Programming: loop versus vector oriented

2010-09-18 Thread Jan private
Hello Petr,

thank you for your ideas. The split() looks most realistic. 

What about this idea:

1. Define three functions Refun1, Refun2, Refun3 for the three different
sections of the calculations (same as you suggested)
2. lambda = (Re = 2320) * Refun1(Re)  + ((Re  2320)  (Re  65 * dk))
* Refun2(Re) etc.

But my thought is that probably the values of RefunXYZ will be
calculated for every value of Re, even if the condition (Re = 2320) is
FALSE (= 0). So that would give a lot of unnecessary function
evaluations.

Regards,
Jan

__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/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] Odp: Programming: loop versus vector oriented

2010-09-17 Thread Jan private
Hello Petr,

 but I think this is how your code really works. Did you try it?

it does, but the R documentation says somewhere:
Warning: for() loops are used in R code much less often than in
compiled languages. Code that takes a `whole object' view is likely to
be both clearer and faster in R.

So I am wondering in what way the whole object view could be applied
to my function.

Best regards,
Jan

  The function should work on the nth elements of the two input vectors
  and put the result into the nth element of the output vector. So it
  would work like c - a + b, only instead of '+' there are more complex
  calculations.

__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/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] Odp: Programming: loop versus vector oriented

2010-09-16 Thread Jan private
Hello Petr,

 If you want to get results of your function for a vector of reynolds and 
 dk you can use function outer and probably get rid of for cycle in the 
 function.
 
 outer(c(100, 530,2410), c(10, 150,200),lambda_wall)
   [,1]   [,2]   [,3]
 [1,] 0.640 0.6400 0.6400
 [2,] 0.1207547 0.12075472 0.12075472
 [3,] 0.1081338 0.04515774 0.04515774

that gives me an array as an answer (and does more calculations than
necessary, in this case 9 instead of 3). The result should be a vector.

The function should work on the nth elements of the two input vectors
and put the result into the nth element of the output vector. So it
would work like c - a + b, only instead of '+' there are more complex
calculations.

Best regards,
Jan

__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/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] Programming: loop versus vector oriented

2010-09-15 Thread Jan private
Dear all,

I am new to R and to it's programming philosophy. The following function
is supposed to work on a vector, but I can't figure out how to do that
without looping through every element of it. Is there a more elegant
way?

Note: I have shortened it, so it is NOT correct from the pipe hydraulics
point of view

# Calculate wall friction factor
# reynolds: Reynolds number
# dk: relative roughness of pipe
lambda_wall - function (reynolds, dk) {
  z - 1
  result - 0

  for (Re in reynolds) {
if (Re = 2320) {
  # Laminar flow
  lambda - 64/Re
} else if (Re  65 * dk[z]) {
  # Turbulent flow
  if (Re  1e+5) {
lambda - 0.3164 / sqrt(sqrt(Re))
  } else {
lambda - 0.309/(log10(Re/7))^2
  }
} else {
  # Intermediate area
  lambdanew - 1 / (2 * log10(3.71 * dk[z]))^2 # Start value
  iter - 0

  repeat {
lambda - lambdanew
lambdanew - 1 / (2 * log10(2.51/(Re * sqrt(lambda)) + 0.27/dk[z]))^2
iter - iter + 1
if ((abs(lambdanew - lambda)  0.001) || (iter  100)) break
  }

  lambda = lambdanew
}

result[z] - lambda
z - z + 1
  }

  result
} # lambda_wall()

__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/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] Calculating with tolerances (error propagation)

2010-09-09 Thread Jan private
Hello Bernardo,

-
If I understood  your problem this script solve your problem:

q-0.15 + c(-.1,0,.1)
h-10 + c(-.1,0,.1)
5*q*h
[1]  2.475  7.500 12.625
-

OK, this solves the simple example. 
But what if the example is not that simple. E.g.

P = 5 * q/h

Here, to get the maximum tolerances for P, we need to divide the maximum
value for q by the minimum value for h, and vice versa. Is there any way
to do this automatically, without thinking about every single step?

There is a thing called interval arithmetic (I saw it as an Octave
package) which would do something like this.

I would have thought that tracking how a (measuring) error propagates
through a complex calculation would be a standard problem of
statistics?? In other words, I am looking for a data type which is a
number with a deviation +- somehow attached to it, with binary operators
that automatically knows how to handle the deviation.

Thank you,  
Jan

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

2010-09-08 Thread Jan private
Dear list,

I am from an engineering background, accustomed to work with tolerances.

For example, I have measured

Q = 0.15 +- 0.01 m^3/s
H = 10 +- 0.1 m

and now I want to calculate

P = 5 * Q * H 

and get a value with a tolerance +-

What is the elegant way of doing this in R?

Thank you,
Jan

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