Plotting Continuous Functions

I see that I've been immortalized on the SciPy MatPlotLib Cookbook web
page for my enquiry on plotting S- and Z-curves. The Boltzman function
serves very well for that purpose, and I've tweaked the example code to
allow me to pass in the two endpoints and the midpoint for each of these
curves.

   Now I need to plot normal curves (a.k.a. Gaussian or bell curves,
depending on the background of the speaker/writer). I see that SciPy has a
class for the normal curve in its stats package, and that the curve shape is
defined by the mean and standard deviation.

   My need is to draw these curves based on the midpoint (== mean) and tail
endpoints (which are not the same as the s.d.). I can copy -- or call -- the
SciPy class from my application code, but I don't know if this is the most
parsimonious approach. Your thoughts are appreciated.

   These plots will be used in two different media, and it two forms. Within
the wxPython application I want to display each curve on one panel, and the
set of curves related to one linguistic variable on a second panel. Then I
need to have the set incorporated into a ReportLab report in .pdf.

Rich

···

--
Richard B. Shepard, Ph.D. | Integrity Credibility
Applied Ecosystem Services, Inc. | Innovation
<http://www.appl-ecosys.com> Voice: 503-667-4517 Fax: 503-667-8863

   Now I need to plot normal curves (a.k.a. Gaussian or bell curves,
depending on the background of the speaker/writer). I see that SciPy has a
class for the normal curve in its stats package, and that the curve shape is
defined by the mean and standard deviation.

For parsimony, I think you're probably best off just using the
Gaussian equation:

def fwhm2k(fwhm):
    '''converts fwhm value to k (see above)'''
    return fwhm/(2 * n.sqrt( n.log( 2 ) ) )

def gauss1d(r, fwhm, c):
    '''returns the 1d gaussian given by fwhm (full-width at half-max),
    and c (centre) at positions given by r
    '''
    return exp( -(r-c)**2 / fwhm2k( fwhm )**2 )

(released to public domain)

   My need is to draw these curves based on the midpoint (== mean) and tail
endpoints (which are not the same as the s.d.).

The midpoint here is c.

It's not clear what you mean by endpoints - if you mean you want to be
able to specify the y value at a given x delta-x away from c, then it
should be relatively simple to solve the equation to find the required
full-width at half-max to achieve these end-points. After a very quick
look (i.e. definitely needs verification), I think

k = sqrt( -(R-c)**2/log(Y) )

where Y is the desired value at distance R-c from the centre.

Your thoughts are appreciated.

I hope that's what you're after.

Angus,

···

On 23/11/2007, Rich Shepard <rshepard@...695...> wrote:
--
AJC McMorland, PhD Student
Physiology, University of Auckland

For parsimony, I think you're probably best off just using the
Gaussian equation:

def fwhm2k(fwhm):
   '''converts fwhm value to k (see above)'''
   return fwhm/(2 * n.sqrt( n.log( 2 ) ) )

def gauss1d(r, fwhm, c):
   '''returns the 1d gaussian given by fwhm (full-width at half-max),
   and c (centre) at positions given by r
   '''
   return exp( -(r-c)**2 / fwhm2k( fwhm )**2 )

   Thank you, Angus. I'll look at the Gaussian explanation to understand the
input values.

The midpoint here is c.

   OK.

It's not clear what you mean by endpoints - if you mean you want to be
able to specify the y value at a given x delta-x away from c, then it
should be relatively simple to solve the equation to find the required
full-width at half-max to achieve these end-points. After a very quick
look (i.e. definitely needs verification), I think

   What I mean is the x value where the tails of the curve have y == 0.0.
These curves are defined by the range of x over which they are valid, and
assume the midpoint is where y == 1.0 (the maximum value). The inflection
points are at y = 0.5; in rare situations that may change.

I hope that's what you're after.

   I'll look at it in detail tomorrow (my time) and the weekend. I, too, hope
that it's what I need.

Much appreciated,

Rich

···

On Fri, 23 Nov 2007, Angus McMorland wrote:

--
Richard B. Shepard, Ph.D. | Integrity Credibility
Applied Ecosystem Services, Inc. | Innovation
<http://www.appl-ecosys.com> Voice: 503-667-4517 Fax: 503-667-8863

Rich Shepard wrote:

···

On Fri, 23 Nov 2007, Angus McMorland wrote:

For parsimony, I think you're probably best off just using the
Gaussian equation:

def fwhm2k(fwhm):
   '''converts fwhm value to k (see above)'''
   return fwhm/(2 * n.sqrt( n.log( 2 ) ) )

def gauss1d(r, fwhm, c):
   '''returns the 1d gaussian given by fwhm (full-width at half-max),
   and c (centre) at positions given by r
   '''
   return exp( -(r-c)**2 / fwhm2k( fwhm )**2 )
    
   Thank you, Angus. I'll look at the Gaussian explanation to understand the
input values.

The midpoint here is c.
    
   OK.

It's not clear what you mean by endpoints - if you mean you want to be
able to specify the y value at a given x delta-x away from c, then it
should be relatively simple to solve the equation to find the required
full-width at half-max to achieve these end-points. After a very quick
look (i.e. definitely needs verification), I think
    
   What I mean is the x value where the tails of the curve have y == 0.0.
These curves are defined by the range of x over which they are valid, and
assume the midpoint is where y == 1.0 (the maximum value). The inflection
points are at y = 0.5; in rare situations that may change.
  
Rich: The tails of a Gaussian never reach zero - they just asymptote to zero for large x.

-Jeff

--
Jeffrey S. Whitaker Phone : (303)497-6313
NOAA/OAR/CDC R/PSD1 FAX : (303)497-6449
325 Broadway Boulder, CO, USA 80305-3328

Jeff,

   For all practical purposes, that's fine. Usually any y value > 0.20 (the
default) is considered functionally equivalent to zero. If the display looks
close enough to the x axis, that fulfills that purpose.

   Because these curves are used to determine degree of set membership the
tails are of no practical value.

Thanks,

Rich

···

On Fri, 23 Nov 2007, Jeff Whitaker wrote:

Rich: The tails of a Gaussian never reach zero - they just asymptote to zero for large x.

--
Richard B. Shepard, Ph.D. | Integrity Credibility
Applied Ecosystem Services, Inc. | Innovation
<http://www.appl-ecosys.com> Voice: 503-667-4517 Fax: 503-667-8863

It's been decades since I've needed to write any statistical or
distributional function code; I used SysStat in DOS for quite some time and
R for the past decade with linux. There's a great difference between using
working code and writing one's one code. :slight_smile:

   Looking more carefully at the SciPy stats module, I see the normal
distribution is included as a class. However, I don't understand all the
functions in that class -- other than normal.pdf(x) -- and do not see where
loc and scale (representing mean and std. deviation) are used. Toward
understanding how to use this code to draw the curves I need, I extracted
the one method into a stand-alone file and tried to plot a Normal/Bell curve
from 0 to 100:

import matplotlib.numerix as nx
import pylab as p
from math import *

def normal(x):
   return exp(-x**2.0/2.0)/sqrt(2.0*pi)

x = nx.arange(0, 100, 0.1)
N = normal(x)
p.plot(x, N, color='red', lw=2)
p.axis([0, 100, 0.0, 1.0])
p.xlabel('Universe of Discourse')
p.ylabel('Membership Grade')
p.show()

   The error returned is:

Traceback (most recent call last):
   File "normal-curve.py", line 17, in ?
     N = normal(x)
   File "normal-curve.py", line 14, in normal
     return exp(-x**2.0/2.0)/sqrt(2.0*pi)
TypeError: only length-1 arrays can be converted to Python scalars

   A clue stick to the meaning of this error message, and how to fix the
problem is needed.

Rich

···

On Thu, 22 Nov 2007, Rich Shepard wrote:

For parsimony, I think you're probably best off just using the
Gaussian equation:

def fwhm2k(fwhm):
   '''converts fwhm value to k (see above)'''
   return fwhm/(2 * n.sqrt( n.log( 2 ) ) )

def gauss1d(r, fwhm, c):
   '''returns the 1d gaussian given by fwhm (full-width at half-max),
   and c (centre) at positions given by r
   '''
   return exp( -(r-c)**2 / fwhm2k( fwhm )**2 )

  Thank you, Angus. I'll look at the Gaussian explanation to understand the
input values.

--
Richard B. Shepard, Ph.D. | Integrity Credibility
Applied Ecosystem Services, Inc. | Innovation
<http://www.appl-ecosys.com> Voice: 503-667-4517 Fax: 503-667-8863