L(x)=c x10e-3553 x
Let's plot it, in S-Plus and also in Mathematica, in several ways, usually with c=1.
First, S-Plus; start it within Emacs or, for command-line die-hards, try "Splus -e" from the command line. Then, once S-Plus is running, type or cut-and-paste:
motif();
help.start();
x <- seq(0,0.01,,101);
f <- x^10 * exp(-3553*x);
plot(x,f,type="l");
f <- function(x) { x^10 * exp(-3553*x); }
plot(x,f(x),type="l",col=3,xlab="Prob",ylab="Density");
title("Example");
text(0.007,.7*max(f(x)), "Lab 3: Exponential Survival Model");
q();
Try to figure out two different built-in S-Plus functions (they will be density functions for different distributions) that lead to the same plot as that of f(), but a different vertical scaling.
Now let's try Mathematica. To begin a mathematica session simply type "math"; on ACPUB, mathematica will respond with its input prompt "In[1]:=". Remember, you need to end each line with a SHIFT RETURN to make Mathematica execute the line.
Plot[x^10 * Exp[-3553*x], {x,0,0.01}]
f[x_] := x^10 * Exp[-3553*x]
ss = FindMinimum[-f[Exp[s]], {s,-5}] (* Now Exp[s] is MLE *)
x = Exp[s] /. ss[[2]] (* Now x is MLE *)
Exp[-365*x] (* Point estimate of p *)
(* SO, p ~ 0.357943 *)
con = Integrate[f[x],{x,0,Infinity}]
N[con]
myplot = Plot[ f[x]/con, {x,0,0.01} ]
Display["lab3.ps",myplot]
Display["!lpr -Psoclp2",myplot]
Integrate[Exp[-x*365]*f[x]/con,{x,0,Infinity}] (* What is this? *)
Quit (* Done! *)