[R] overlay lattice histograms with goodness-of-fit pdfs

[Brad Christoffersen]

Mange tak!

FYI, this is the way it is able to run (I was going to attach
"station.precip.R", but I read that attaching files is not recommended 
- let me
know if you would like it)

x <- dget(file="C://Documents and Settings/Bradley/My
Documents/Arizona/CourseResources/ATMO529/station.precip.R")
histogram(~ data | month * station, data = sta.stack,
	subset = type=="precip" & month %in% c("Dec","Jan","Feb"),
	xlab = "Precipitation (mm)",
	type = "density",
	panel = function(x, ...) {
		panel.histogram(x, ...)
		panel.mathdensity(dmath = dnorm, col = "black",
					args = list(mean = mean(sta.stack$data), sd = sd(sta.stack$data)))
		panel.mathdensity(dmath = dgamma, col = "black",
					args = list(shape = (mean(sta.stack$data))^2 / (stdev(sta.stack$data))^2,
							scale = (stdev(sta.stack$data))^2 / mean(sta.stack$data)))
	})

Now, what would be great is to be able to reference the different calls to
panel.mathdensity() so that the corresponding probability * histogram area ( =
counts) in each bin can be used to compute a simple chi-square 
goodness-of-fit.
  I tried calling panel.mathdensity() outside of histogram(), but I don't think
this is right - it returns NULL.  I also looked at chisq.test, but this 
doesn't
support trellis formulas.  Any thoughts or leads?

Thanks,
Brad Christoffersen



Quoting Frede Aakmann Tøgersen <FredeA.Togersen at agrsci.dk>:

> The following is one of the examples in the help page for histogram:
>
>      histogram( ~ height | voice.part, data = singer,
>                xlab = "Height (inches)", type = "density",
>                panel = function(x, ...) {
>                    panel.histogram(x, ...)
>                    panel.mathdensity(dmath = dnorm, col = "black",
>                                      args = list(mean=mean(x),sd=sd(x)))
>                } )
>
> This should give you some thing to start from.
>
> Also using the subset argument of the lattice functions will make 
> make your code more readable. Instead of your code
>
> histogram(~ data | month * station,
> 	data = sta.stack[sta.stack[,"type"]=="precip" & 
> (sta.stack[,"month"]=="Dec" | sta.stack[,"month"]=="Jan" | 
> sta.stack[,"month"]=="Feb"),],
> 	xlab = "Precipitation (mm)")
>
> you can use (not tested because you didn't supply a reproducable example)
>
> histogram(~ data | month * station, data = sta.stack
> 	subset = type==precip & month %in% c("Dec", "Jan", "Feb"),
> 	xlab = "Precipitation (mm)")
>
>
> Med venlig hilsen
> Frede Aakmann Tøgersen
>
>
>
>
>> -----Oprindelig meddelelse-----
>> Fra: r-help-bounces at stat.math.ethz.ch
>> [mailto:r-help-bounces at stat.math.ethz.ch] På vegne af Brad
>> Christoffersen
>> Sendt: 10. september 2007 12:08
>> Til: R-help at stat.math.ethz.ch
>> Emne: [R] overlay lattice histograms with goodness-of-fit pdfs
>>
>> Hello,
>>
>> I am new to R exploratory data analysis and plotting.  Is
>> anyone aware of a way to overlay a set of conditional
>> histograms with conditional PDFs?  Below, I generate a
>> lattice plot of precipitation histograms based on different
>> months and stations, given a subset of the dataset:
>>
>>
>> histogram(~ data | month * station,
>> 	data = sta.stack[sta.stack[,"type"]=="precip" &
>> (sta.stack[,"month"]=="Dec" | sta.stack[,"month"]=="Jan" |
>> sta.stack[,"month"]=="Feb"),],
>> 	xlab = "Precipitation (mm)")
>>
>>
>> I previously used a combination of the low-level 'lines()'
>> and 'dgamma()'
>> functions to overlay a gamma PDF onto a single histogram.
>> Now what I would like to do is to do the same thing, but with
>> a function that allows me to specify a formula similar to
>> that in the histogram function above
>>
>> [SomeKindOfPDF] ~ [x-range] | month * station
>>
>> which will plot the PDF with the appropriate factors (month
>> and station).
>>
>> All I'm looking for is for someone to get me going in the
>> right direction with a useful package or function to use.
>>
>> Any help is much appreciated!
>> Brad Christoffersen
>>
>> ______________________________________________
>> R-help at stat.math.ethz.ch 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.
>>