742.sagews

# In this notebook we create a bimodal pdf by concatenating together two normal pdfs (see http://www.talkstats.com/showthread.php?t=9776)

# In these following steps the constants necessary so that the function is a pdf are found

n(integral( exp( -(x+1)^2  / 2 ), x, -Infinity, 0))

2.10893852920804

k1 = (1/2) * 1 / n(integral( exp(-(x+1)^2 / 2 ), x, -Infinity, 0)); k1

0.237086094770037

n(integral( exp( -(x-1)^2 / 2), x, 0, Infinity))

2.10893852920804

k2 = (1/2) * 1 / n(integral( exp(-(x-1)^2 / 2 ), x, 0, Infinity)); k2

0.237086094770037

# Since it is symmetric k1 and k2 are the same

k = k1

n(k * n(integral(exp( -(x+1)^2 / 2 ), x, -Infinity, 0)) + \
  k * n(integral(exp( -(x-1)^2 / 2 ), x, 0, Infinity )) )

1.00000000000000

def bimodal_pdf(x):
    if ( x > -Infinity and x < 0 ):
        return k * exp( -(x+1)**2 / 2)
    else:
        return k * exp( -(x-1)**2 / 2)

# note: this method of doing a piecewise function seems to be better suited for plotting than using the piecewise code found in the sage documentation.

plot( bimodal_pdf, -5, 5)

# The modes for this bimodal pdf are -1 and 1 (the mode being defined for a continuous pdf as the value which maximizes the pdf function)