[R] looping problem

Roslina Zakaria zroslina at yahoo.com
Mon Apr 14 05:56:41 CEST 2008


Hi R-users,

I would like to do looping for this process below to estimate alpha beta from gamma distribution:
Here are my data:
day_data1 <- 
        1    2    3    4    5    6     7    8    9   10   11   12   13   14  15   16   17   18   19   20   21   22   23
1943 48.3 18.5  0.0  0.0 18.3  0.0   0.0  0.0  0.0  0.0  0.0  0.0  2.0  0.0 0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.8  2.8
1944  0.0  0.0  0.0  0.0  0.0  0.0   0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0 0.0  0.0  6.6  0.0  0.0  0.0  0.0  0.0  0.0
1945  5.3  0.0  0.0  0.0  0.0  2.5   0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0 0.0  0.0  5.3  0.0  0.0  0.0  0.0  0.0  0.0
1946  0.0  0.0  0.0  0.0  0.0  2.3   0.0  0.0  0.0  0.0  4.8  0.3  1.5  0.0 0.8  0.0  0.0  5.8 70.6 12.4  0.5 23.6  0.0
1947  0.0  0.0  0.0  0.0  0.0  0.0   0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0 0.0  2.0  0.0  0.0  0.0  2.8  0.0  0.0  0.0
1948  0.3 20.1  0.0  0.0  0.0  0.0   0.0  0.0  0.0  0.0  1.5  0.5  0.0  0.0 0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0

## To extract all the positive values
x1 <- day_data1[,1]  
x2 <- day_data1[,2]
x3 <- day_data1[,3]
x4 <- day_data1[,4]
x5 <- day_data1[,5]
tol <- 1E-6  
a1 <- x1[x1>tol]
a2 <- x2[x2>tol]
a3 <- x3[x3>tol]
a4 <- x4[x4>tol]
a5 <- x5[x5>tol]

library(MASS)
## Example January Charleville 1943-2007
fitdistr(a1,dgamma, list(shape = 1, rate = 0.1), lower = 0.01)
fitdistr(a2,dgamma, list(shape = 1, rate = 0.1), lower = 0.01)
fitdistr(a3,dgamma, list(shape = 1, rate = 0.1), lower = 0.01)
fitdistr(a4,dgamma, list(shape = 1, rate = 0.1), lower = 0.01)
fitdistr(a5,dgamma, list(shape = 1, rate = 0.1), lower = 0.01)

Here is my code: 

alpha.beta <- function(data,n)
{  tol <- 1E-6
  { for (i in 1:n)
    xi <- data[,i]
    ai <- xi [xi > tol]
    fit <- fitdistr(ai,dgamma, list(shape = 1, rate = 0.1), lower = 0.01)
  }
  fit
}

I’m not sure what went wrong since it gives only one output by right 31 outputs.
Thank you for your attention. 

__________________________________
t spam protection around 
http://mail.yahoo.com



More information about the R-help mailing list