Post date: Jul 05, 2021 12:31:26 PM
R script to find example variance-extremal distributions for four configurations of Dempster-Shafer structures having only two interval focal elements
####################################################################################
# population statistics are computed by dividing by n, not n-1
pvar = function(x) (1 - 1/length(x)) * var(x)
psd = function(x) sqrt(1 - 1/length(x)) * sd(x)
specific = function(a) runif(1) * (right(a) - left(a)) + left(a)
endpoint = function(a) if (runif(1) < 0.5) return(left(a)) else return(right(a))
cfg = function(a, how=endpoint) {
A = NULL
for (i in 1:length(a)) A = c(A, specific(a[i]))
A
}
search = function(a, m=Inf, M=-Inf, smallest=NULL, largest=NULL, some=1, many=10000, how=endpoint) {
A = NULL
for (k in 1:many) {
for (n in 1:some) A = c(A,cfg(a,how))
pv = pvar(A)
if (pv < m) {m = pv; smallest = A}
if (M < pv) {M = pv; largest = A}
}
#cat('m,M',m,' ',M,'\n')
list(m=m,M=M,smallest=smallest,largest=largest)
}
finetune = function(a,m,extremal,f=`<`,dispersion=0.05,many=1000) {
n = length(a)
r = rep(0,n)
for (i in 1:n) r[[i]] = (extremal[[i]] - a[[i]]@lo) / (a[[i]]@hi - a[[i]]@lo)
for (k in 1:many) {
r = rnorm(n,r,dispersion)
r = ifelse(r<0, abs(r), r)
r = ifelse(1<r, 2-r, r)
if (any(r<0 | 1<r)) stop('Too dispersed')
A = NULL
for (i in 1:n) A = c(A, r[[i]] * (a[[i]]@hi - a[[i]]@lo) + a[[i]]@lo)
pv = pvar(A)
if (f(pv,m)) {m = pv; extremal = A}
}
if (m<1e-5) m = 0 ##################################
list(m=m,extremal=extremal)
}
explore = function(a, many=1000) {
m = Inf
M = -Inf
smallest=NULL
largest=NULL
A = c(a[[1]]@lo,a[[2]]@lo); pv=pvar(A); if (pv<m) {m=pv; smallest=A}; if (M<pv) {M=pv; largest=A}
A = c(a[[1]]@lo,a[[2]]@hi); pv=pvar(A); if (pv<m) {m=pv; smallest=A}; if (M<pv) {M=pv; largest=A}
A = c(a[[1]]@hi,a[[2]]@lo); pv=pvar(A); if (pv<m) {m=pv; smallest=A}; if (M<pv) {M=pv; largest=A}
A = c(a[[1]]@hi,a[[2]]@hi); pv=pvar(A); if (pv<m) {m=pv; smallest=A}; if (M<pv) {M=pv; largest=A}
A = c(a[[1]]@lo,a[[2]]@lo,a[[1]]@hi,a[[2]]@hi); pv=pvar(A); if (pv<m) {m=pv; smallest=A}; if (M<pv) {M=pv; largest=A}
for (s in 1:4) for (f in c(endpoint, specific)) {
e = search(a,m,M,smallest,largest,s,many,f)
m = e$m
M = e$M
smallest = e$smallest
largest = e$largest
}
e = finetune(a,m,smallest,f=`<`,dispersion=0.05,many=2000)
m = e$m
smallest = e$extremal
e = finetune(a,M,largest,f=`>`,dispersion=0.05,many=2000)
M = e$m
largest = e$extremal
n = length(a)
c = NULL
for (i in 1:n) c = c(c, a[[i]]@lo,a[[i]]@hi)
plot(NULL, xlim=range(c),ylim=c(0,n+1))
for (i in 1:n) lines(c(a[[i]]@lo,a[[i]]@hi),c(i,i),lwd=3)
points(smallest,rep(1:n,length(smallest)/n),col='red',cex=2,lwd=3)
points(largest,rep(1:n,length(largest)/n),col='blue',cex=3,,lwd=3)
title(paste('variance = [',signif(m,2),',',signif(M,2),']'))
list(v = c(m,M), smallest = smallest, largest = largest)
}
####################################################################################
par(mfrow=c(1,1))
xx = seq(5,6,length.out=100)
four = two = NULL
for (x in xx) {
four = c(four, pvar(c(1,2,4,x)))
two = c(two, pvar(c(2,x)))
}
plot(xx,four,col='red',lwd=2)
points(xx,two,col='blue',lwd=2)
####################################################################################
par(mfrow=c(2,2))
explore(c(interval(1, 2), interval(3, 4)))
explore(c(interval(1, 3), interval(2, 4)))
explore(c(interval(1, 5), interval(2, 4)))
explore(c(interval(1, 6), interval(2, 4)))
par(mfrow=c(1,1))
####################################################################################