splice <- function(c, fns, range=c(0,Inf)) {
n = length(fns)
breaks = c(range[1], c, range[2])
pcw <- function(x) {
g<-cut(x, breaks)
unsplit(Map(function(f,x) f(x), fns, split(x, g)), g)
}
integrals = as.numeric(lapply(1:n, function(i) {force(i); integrate(pcw, breaks[i], breaks[i+1])$value }))
ratios = as.numeric(lapply((1:(n-1)), function(i) {force(i); fns[[i+1]](c[i])/fns[[i]](c[i]) }))
a = numeric(n)
products = as.numeric(lapply( (1:(n-1)), function(i) {force(i); prod(ratios[i:(n-1)])}))
a[n] = (sum(products*integrals[1:length(integrals)-1]) + integrals[length(integrals)])^(-1)
for(i in (n-1):1) a[i] = a[i+1]*ratios[i] # backwards recursion
fnsWeighted = lapply(1:n, function(i) {force(i); function(x) f[[i]](x)*a[i]})
function(x) {
g<-cut(x, breaks)
unsplit(Map(function(f,x) f(x), fnsWeighted, split(x, g)), g)
}
}
splice.plot <- function(fns, g, c, range=c(0, 100), title="lnorm(1, 0.9078), exp(1)", ylim=c(0,0.1)) {
n = length(fns)
breaks = c(range[1], c, range[2])
plot(g, range[1], range[2], main=title, ylim=ylim)
lapply(1:n, function(i) {force(i); plot(fns[[i]], breaks[i], breaks[i+1], col=i+1, lty=3, add=TRUE) })
lapply((1:(n-1)), function(i) {force(i); plot(fns[[i]], breaks[i+1], range[2], col=i+1, add=TRUE)})
lapply(n:2, function(i) {force(i); plot(fns[[i]], range[1], breaks[i], col=i+1, add=TRUE) })
abline(v=c, lty=3)
}
x <- readline("2 pdfs. Press any key to continue")
## 2 pdfs. Press any key to continue
c = 2
f = c(function(x) dlnorm(x, meanlog=1, sdlog=0.9078),
function(x) dexp(x, rate=1))
g = splice(c, f) # spliced pdf
splice.plot(f, g, c, range=c(0, 10), title="lnorm(1, 0.9078), exp(1)", ylim=c(0, 0.5))

x <- readline("3 pdfs. Press any key to continue")
## 3 pdfs. Press any key to continue
c = c(20, 40)
f = c(function(x) dlnorm(x, meanlog=3.5, sdlog=0.9078),
function(x) dweibull(x, shape=3, scale=100/sqrt(pi)),
function(x) dgamma(x, shape=50, scale=1))
g = splice(c, f)
splice.plot(f, g, c, range=c(0, 100), title="lnorm(3.5, 0.9078), weibull(3, 56.42), gamma(50, 1)", ylim=c(0, 0.06))

x <- readline("4 pdfs. Press any key to continue")
## 4 pdfs. Press any key to continue
c = c(20, 40, 60)
f = c(function(x) dlnorm(x, meanlog=3.5, sdlog=0.9087),
function(x) dweibull(x, shape=3, scale=100/sqrt(pi)),
function(x) dgamma(x, shape=50, scale=1),
function(x) dweibull(x, shape=5, scale=75))
g = splice(c, f)
splice.plot(f, g, c, range=c(0, 100), title="lnorm(3.5, 0.9078), weibull(3, 56.42), gamma(50, 1), weibull(5,75)", ylim=c(0, 0.06))
