- fun = function(x,mu,sigma,alpha, beta){
- re = function(y){
- return(exp(-(y-mu)^2/(2*sigma^2))*exp(alpha+beta*y)/(1+exp(alpha+beta*y)))
- }
- return(re(x)/(integrate(re, -100, 100))$value)
- }
- funv = Vectorize(fun)
- x = seq(-10, 10, len=200)
- y1=funv(x,mu=0,sigma=1,alpha=1,beta=1)
- y2=funv(x,mu=0,sigma=1,alpha=1,beta=-1)
- plot(x, y1, type = 'n', xlab = 'x', ylab = 'y')
- lines(x, y1, lty=2, lwd=2)
- lines(x, y2, lty=3, lwd=2)
- legend('topright', c('y1','y2'), lty=2:3, lwd=2, bty='n')
复制代码
|