[R] Why two curves and numerical integration look so different?
C W
tmrsg11 at gmail.com
Fri Feb 12 17:29:16 CET 2016
Hi Peter,
Great, let me try that and get back to you on my findings in a few hours!
:)
On Fri, Feb 12, 2016 at 11:09 AM, peter dalgaard <pdalgd at gmail.com> wrote:
> I don't see here FAQ 7.31 comes in either (for once!...)
>
> However, either the density is unnormalized and the integral is not 1, or
> the integral is 1 and it is normalized. The one in the picture clearly does
> not integrate to one. You can fit a rectangle of size 0.1 by 1e191 under
> the curve so the integral should be > 1e190 .
>
> As depicted, I don't see why a plain integral from .5 to 1.5 shouldn't
> work.
>
> -pd
>
> On 12 Feb 2016, at 16:57 , C W <tmrsg11 at gmail.com> wrote:
>
> > Hi Bert,
> >
> > Yay fantasyland!
> >
> > In all seriousness, You are referring to this?
> >
> https://cran.r-project.org/doc/FAQ/R-FAQ.html#Why-doesn_0027t-R-think-these-numbers-are-equal_003f
> >
> > In particular, you mean this: .Machine$double.eps ^ 0.5?
> >
> > Thanks!
> >
> > On Fri, Feb 12, 2016 at 10:53 AM, Bert Gunter <bgunter.4567 at gmail.com>
> > wrote:
> >
> >> You are working in fantasyland. Your density is nonsense.
> >>
> >> Please see FAQ 7.31 for links to computer precision of numeric
> >> calculations.
> >>
> >>
> >> Cheers,
> >> Bert
> >> Bert Gunter
> >>
> >> "The trouble with having an open mind is that people keep coming along
> >> and sticking things into it."
> >> -- Opus (aka Berkeley Breathed in his "Bloom County" comic strip )
> >>
> >>
> >> On Fri, Feb 12, 2016 at 7:36 AM, C W <tmrsg11 at gmail.com> wrote:
> >>> Hi David,
> >>>
> >>> This is the Gaussian looking distribution I am trying to integrate.
> >>>
> >>
> https://drive.google.com/file/d/0B2xN0-A6iTB4NThIZ2tYdGxHc00/view?usp=sharing
> >>>
> >>> Notice the unnormalized density goes up as high as 2.5*101^191.
> >>>
> >>> I tried to create small intervals like
> >>>> seq(0.5, 1.3, by = 10^(-8))
> >>>
> >>> but that doesn't seem to be good enough, as we know, it should
> integrate
> >> to
> >>> 1.
> >>>
> >>> On Thu, Feb 11, 2016 at 3:32 PM, David Winsemius <
> dwinsemius at comcast.net
> >>>
> >>> wrote:
> >>>
> >>>>
> >>>>> On Feb 11, 2016, at 11:30 AM, C W <tmrsg11 at gmail.com> wrote:
> >>>>>
> >>>>> Hi David,
> >>>>>
> >>>>> My real function is actually a multivariate normal, the simple toy
> 1-d
> >>>> normal won't work.
> >>>>>
> >>>>> But, you gave me an idea about restricting the bounds, and focus
> >>>> integrating on that. I will get back to you if I need any further
> >>>> assistance.
> >>>>
> >>>> You'll probably need to restrict your bounds even more severely than I
> >> did
> >>>> in the 1-d case (using 10 SD's on either side of the mean) . You might
> >> need
> >>>> adequate representation of points near the center of your
> >> hyper-rectangles.
> >>>> At least that's my armchair notion since I expect the densities tail
> off
> >>>> rapidly in the corners. You can shoehorn multivariate integration
> around
> >>>> the `integrate` function but it's messy and inefficient. There are
> other
> >>>> packages that would be better choices. There's an entire section on
> >>>> numerical differentiation and integration in CRAN Task View: Numerical
> >>>> Mathematics.
> >>>>
> >>>> --
> >>>> David.
> >>>>
> >>>>
> >>>>>
> >>>>> Thank you so much!
> >>>>>
> >>>>> On Thu, Feb 11, 2016 at 2:06 PM, David Winsemius <
> >> dwinsemius at comcast.net>
> >>>> wrote:
> >>>>>
> >>>>>> On Feb 11, 2016, at 9:20 AM, C W <tmrsg11 at gmail.com> wrote:
> >>>>>>
> >>>>>> I want to do numerical integration w.r.t. mu: P(mu) × N(mu, 0.00001)
> >>>>>>
> >>>>>> Because the variance is small, it results in density like:
> >> 7.978846e+94
> >>>>>>
> >>>>>> Is there any good suggestion for this?
> >>>>>
> >>>>> So what's the difficulty? It's rather like the Dirac function.
> >>>>>
> >>>>>> integrate( function(x) dnorm(x, sd=0.00001), -.0001,0.0001)
> >>>>> 1 with absolute error < 7.4e-05
> >>>>>
> >>>>>
> >>>>> --
> >>>>> David.
> >>>>>
> >>>>>>
> >>>>>> Thanks so much!
> >>>>>>
> >>>>>>
> >>>>>> On Thu, Feb 11, 2016 at 9:14 AM, C W <tmrsg11 at gmail.com> wrote:
> >>>>>>
> >>>>>>> Wow, thank you, that was very clear. Let me give it some more runs
> >>>> and
> >>>>>>> investigate this.
> >>>>>>>
> >>>>>>> On Thu, Feb 11, 2016 at 12:31 AM, William Dunlap <
> >> wdunlap at tibco.com>
> >>>>>>> wrote:
> >>>>>>>
> >>>>>>>> Most of the mass of that distribution is within 3e-100 of 2.
> >>>>>>>> You have to be pretty lucky to have a point in sequence
> >>>>>>>> land there. (You will get at most one point there because
> >>>>>>>> the difference between 2 and its nearest neightbors is on
> >>>>>>>> the order of 1e-16.)
> >>>>>>>>
> >>>>>>>> seq(-2,4,len=101), as used by default in curve, does include 2
> >>>>>>>> but seq(-3,4,len=101) and seq(-2,4,len=100) do not so
> >>>>>>>> curve(..., -3, 4, 101) and curve(..., -2, 4, 100) will not show
> >> the
> >>>> bump.
> >>>>>>>> The same principal holds for numerical integration.
> >>>>>>>>
> >>>>>>>>
> >>>>>>>> Bill Dunlap
> >>>>>>>> TIBCO Software
> >>>>>>>> wdunlap tibco.com
> >>>>>>>>
> >>>>>>>> On Wed, Feb 10, 2016 at 6:37 PM, C W <tmrsg11 at gmail.com> wrote:
> >>>>>>>>
> >>>>>>>>> Dear R,
> >>>>>>>>>
> >>>>>>>>> I am graphing the following normal density curve. Why does it
> >> look
> >>>> so
> >>>>>>>>> different?
> >>>>>>>>>
> >>>>>>>>> # the curves
> >>>>>>>>> x <- seq(-2, 4, by=0.00001)
> >>>>>>>>> curve(dnorm(x, 2, 10^(-100)), -4, 4) #right answer
> >>>>>>>>> curve(dnorm(x, 2, 10^(-100)), -3, 4) #changed -4 to -3, I get
> >> wrong
> >>>>>>>>> answer
> >>>>>>>>>
> >>>>>>>>> Why the second curve is flat? I just changed it from -4 to -3.
> >>>> There is
> >>>>>>>>> no density in that region.
> >>>>>>>>>
> >>>>>>>>>
> >>>>>>>>> Also, I am doing numerical integration. Why are they so
> >> different?
> >>>>>>>>>
> >>>>>>>>>> x <- seq(-2, 4, by=0.00001)
> >>>>>>>>>> sum(x*dnorm(x, 2, 10^(-100)))*0.00001
> >>>>>>>>> [1] 7.978846e+94
> >>>>>>>>>> x <- seq(-1, 4, by=0.00001) #changed -2 to -1
> >>>>>>>>>> sum(x*dnorm(x, 2, 10^(-100)))*0.00001
> >>>>>>>>> [1] 0
> >>>>>>>>>
> >>>>>>>>> What is going here? What a I doing wrong?
> >>>>>>>>>
> >>>>>>>>> Thanks so much!
> >>>>>>>>>
> >>>>>>>>> [[alternative HTML version deleted]]
> >>>>>>>>>
> >>>>>>>>> ______________________________________________
> >>>>>>>>> R-help at 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.
> >>>>>>>>>
> >>>>>>>>
> >>>>>>>>
> >>>>>>>
> >>>>>>
> >>>>>> [[alternative HTML version deleted]]
> >>>>>>
> >>>>>> ______________________________________________
> >>>>>> R-help at 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.
> >>>>>
> >>>>> David Winsemius
> >>>>> Alameda, CA, USA
> >>>>>
> >>>>>
> >>>>
> >>>> David Winsemius
> >>>> Alameda, CA, USA
> >>>>
> >>>>
> >>>
> >>> [[alternative HTML version deleted]]
> >>>
> >>> ______________________________________________
> >>> R-help at 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.
> >>
> >
> > [[alternative HTML version deleted]]
> >
> > ______________________________________________
> > R-help at 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.
>
> --
> Peter Dalgaard, Professor,
> Center for Statistics, Copenhagen Business School
> Solbjerg Plads 3, 2000 Frederiksberg, Denmark
> Phone: (+45)38153501
> Office: A 4.23
> Email: pd.mes at cbs.dk Priv: PDalgd at gmail.com
>
>
[[alternative HTML version deleted]]
More information about the R-help
mailing list