Take the 2-minute tour ×
Cross Validated is a question and answer site for people interested in statistics, machine learning, data analysis, data mining, and data visualization. It's 100% free, no registration required.

The probability density function I have is as following:

enter image description here

Now, I want to rewrite the cumulative probability function bellow into R.

enter image description here

The variables ($S$, $\mu$, $\sigma$) are constant variables. For $t$ i.e. t <- seq(0.1,1,0.1)

How the cumulative probability would look like in R? Non of my solutions yields desired output!

EDIT:

f1 = function(t,S,sigma,mu){((mu*t+S)/(sigma*sqrt(t)))}
f2 = function(t,S,sigma,mu){((mu*t-S)/(sigma*sqrt(t)))}

1-(pnorm(f1(t,S,sigma,mu),mu,sigma,1,0)-exp(-(2*S*mu/sigma^2))*pnorm(f2(t,S,sigma,mu),mu,sigma,1,0))

This kind of gives me reasonable results but not sure if correct.

share|improve this question
1  
This is not a CDF as it does not converge to 1 if t goes to infinity. –  Michael Mayer Jan 12 at 19:10
    
You are right! What I need is cumulative probability and not CDF! –  Max Jan 12 at 19:11
1  
Max - (i) what do you mean by cumulative probability if not CDF? (ii) Where do your equations come from? Are you sure you have them right? –  Glen_b Jan 12 at 21:00

2 Answers 2

up vote 1 down vote accepted

This would implement your function in R. What part of the output didn't make sense when you tried?

 cumulative = function(t){
      S = somenumber
      mu = somenumber
      sigma = somenumber
      cd = pnorm((mu*t+S)/(sigma*sqrt(t)))-exp(-2*S*mu/sigma^2)*pnorm((mu*t-S)/(sigma*sqrt(t)))
      return(cd)
    }

That being said, what are your parameter values? For combinations I've tried, this seems like a decreasing function of t.

share|improve this answer
    
I'm just trying to replicate someone elses work but don't have the complete workouts so I'm second guessing some things here. –  Max Jan 12 at 20:03
    
My solution (see EDIT) provides me with the desired outcome (same as yours just with the 1-(pnorm...) so something must be missing here. In general your answer is correct. Thanks. –  Max Jan 12 at 20:11
    
@Max Your problem is in the pnorm call. You call pnorm(f1(...),mu,sigma,0,1) so you're not computing from a N(0,1) distribution, but from a N(mu,sigma) distribution. Also by including the 1,0 in your call, you are reversing the calculation, pnorm(2,mu,sigma,0,0)=1-pnorm(2,mu,sigma,1,0) –  ltronneberg Jan 12 at 20:45
    
@Max And we dont get the same results, atleast not when I tried to compare the two. –  ltronneberg Jan 12 at 20:55

A CDF is defined as non-decreasing, right continuous function F, with

  • $\lim_{x\to -\infty}F(x)=0,$ and

  • $\lim_{x\to +\infty}F(x)=1.$

I am just stating the obvious here, but when we let $t<0$ we have already a problem as your real-valued function contains square roots of t. Similarly, we need $\sigma >0$ and perhaps a lot of other things.

In either case, I will just go ahead and assume you know what you are doing, and show you the code you asked for.

# custom cdf
custom.cdf <- function(S,sigma,mu,t)
{
    q1 <- pnorm((mu*t + S)/(sigma*sqrt(t)))
    q2 <- pnorm((mu*t - S)/(sigma*sqrt(t)))
    e.term <- exp(-((2*S*mu)/(sigma^2)))
    res <- q1 - e.term*q2
    return(res)
}    

# set parameters
S <- 1
sigma <- 1
mu <- 2    

t <- seq(0.1,1,0.1)
plot(t,custom.cdf(S,sigma,mu,t),typ="l")    

I have chosen a larger support to see a bit more of the curve. Check the formula (mine and yours) for any mistakes.

t <- seq(1e-6,1,1e-6)
plot(t,custom.cdf(S,sigma,mu,t),typ="l")
share|improve this answer
1  
Shouldn't q1 <- pnorm((mu*t + S)/sigma*sqrt(t)) be pnorm((mu*t+S)/(sigma*sqrt(t)))? As sqrt(t) is also in the denominator –  ltronneberg Jan 12 at 20:03
    
Yes, spot on. I was getting even more spectacular outputs here :) –  Max Jan 12 at 20:05
    
Yes that's correct. Thanks! I will edit and get rid of the graphic. –  means-to-meaning Jan 12 at 20:16

Your Answer

 
discard

By posting your answer, you agree to the privacy policy and terms of service.

Not the answer you're looking for? Browse other questions tagged or ask your own question.