On Jul 16, 2009, at 5:53 PM, KARSTEN G HOLMQUIST wrote:

Hello,

I have a series of questions that I hope will be simple to answer. Basically I would like a code to do the following so that I can compute the distribution free test for the slope of a postulated regression line (Theil test). As I am testing the null hypothesis that slope = 0 against the general alternative the slope does not equal 0, it should be pretty straight-forward.

I have a data frame as follows;


index   Score   Mean_Temp       SD_Mean_Temp
i       10      3.7             0.65
i+1     12      7.5             2.11
i+2     6       10.3            0.85
...     ...     ...             ...
n       8       12.6            4.2

The first question is how can I compute the n(n-1)/2 differences between Score i and Score j where i is from 1 to (n-1) and j is from (i+1) to n? I would then like R to spit out the differences in a file something like the following;

i       j       Difference
1       2       2
1       3       -4
...     ..      ...
(n-1)   n       D(n-1,n)

Without an example, I am going to suggest something:

dfg<-data.frame(index=1:10, Score=sample(1:20, 10),Mean_Temp=10*abs(rnorm(10)), SD_Mean_Temp=abs(rnorm(10)) )

The trick is to get the differences and here is one method:

mapply(FUN="-", dfg$Mean_Temp[combn(dfg$index,2)[1,] ], dfg $Mean_Temp[combn(dfg$index,2)[2,] ])

Should be pretty obvious how to line up the indices since they are just the arguments to "[".



The second question is for each indexed sample from i to n, I would like to use something like rnorm(n, mean = 0, sd = 1) so that I could generate 1000 random draws from the distribution specified so that the arguments for rnorm () are as follows; mean = the value of Mean_Temp & sd = SD_Mean_Temp for indexed values i to n. I would, of course, then like R to spit out the massive table with 1000 columns of randomly generated temperature values for each index i to n.

You can get 1000 such values with rnorm(1000, Mean_Temp, SD_Mean_Temp)

Perhaps something along these lines:

Mtx <- matrix( , ncol=1000, nrow=nrow(dfg))
Mtx[1:nrow(dfg), ] <- apply(dfg, 1, function(x) rnorm(1000, x[3], x[4]) )


Finally I would like to compute the slope estimator associated with the Theil statistic. All I need to do is compute the n(n-1)/n
^^^^^^^^^^==n, is that what you meant? Or is it n(n-1)/2?

individual slope values; Sij = (Yj -Yi)/(Xj-Xi) where the Ys are the Scores for samples i to n and the Xs are the set of randomly generated temperature values, as before, 1 less than or equal to i less than j less than or equal to n. The median is then the estimator of Beta (the slope). I would like to compute the slopes for all 1000 sets of randomly generated temperature values.

There is a package, mblm, that computes the Theil estimator. I don't see why we should reinvent that wheel.

(Found via tsearching on theil statistic)
http://search.r-project.org/cgi-bin/namazu.cgi?query=theil+statistic&max=100&result=normal&sort=score&idxname=functions&idxname=Rhelp08&idxname=views
...and following what seemed like promising links


library(mblm)

--

David Winsemius, MD
Heritage Laboratories
West Hartford, CT

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

Reply via email to