1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194

/* -*- mode: java; c-basic-offset: 4; indent-tabs-mode: nil -*- */

package net.willware.semweb;

import java.util.ArrayList;
import java.util.List;

import com.hp.hpl.jena.graph.Node;
import com.hp.hpl.jena.query.*;
import com.hp.hpl.jena.rdf.model.*;
import com.hp.hpl.jena.reasoner.*;
import com.hp.hpl.jena.reasoner.rulesys.*;

/**
 * Let A and B be random boolean variables, with A being an
 * unobservable cause and B being an observable effect, and Pr(B|A) and
 * Pr(B|~A) are given. This might be the case if the situation has some
 * structure that dictates the conditional probabilities, like a 6.432
 * problem, or it just might be empirical. From these we can compute
 * probabilities for the four cases.
 * <ul><li>P11 = Pr(B^A) = Pr(A) Pr(B|A)</li>
 * <li>P10 = Pr(~B^A) = Pr(A) (1-Pr(B|A))</li>
 * <li>P01 = Pr(B^~A) = (1-Pr(A)) Pr(B|~A)</li>
 * <li>P00 = Pr(~B^~A) = (1-Pr(A)) (1-Pr(B|~A))</li></ul>
 * 
 * Treat (Pr(B|A), Pr(B|~A)) as one piece of information, and Pr(A) as a
 * separate independent piece of information. The first piece reflects
 * your beliefs about the machinery connecting A to B, and the second
 * reflects your beliefs about the likelihood of A, so these are
 * legitimately separate concerns.<p>
 * 
 * If we observe that B=1, then we want to replace Pr(A) with our
 * previous estimate for Pr(A|B), which is given by our old numbers as
 * P11/Pr(A) = P11/(P11+P01), and this becomes our posterior probability
 * for A.<p>
 * 
 * The RDF graph for doing all this looks like this.<p>
 * 
 * <img src="BayesRdfHack.png"><p>
 * 
 * <ul><li><a href="http://willware.blogspot.com/search/label/bayesian">Blog posts</a></li>
 * <li><a href="http://code.google.com/p/wware-autosci/source/browse/bayesnets/updateMath.py">
 *   Bayesian update math</a></li></ul>
 */
public class BayesianInference extends JenaUtil {
    
    private static final String baseUri = "http://127.0.0.1/bayes.rdf#";
    private Model model;
    private Property valueProperty;
    private Property probProperty;

    public static void main(String[] args) {
        BayesianInference bi = new BayesianInference();
        System.out.println("Initial Pr(A) = " + bi.getProbability("A"));
        bi.setBooleanValue("B", true);
        for (int i = 0; i < 5; i++) {
            bi.updateProb();
            System.out.println("T -> Pr(A) = " + bi.getProbability("A"));
        }
        bi.setBooleanValue("B", false);
        for (int i = 0; i < 10; i++) {
            bi.updateProb();
            System.out.println("F -> Pr(A) = " + bi.getProbability("A"));
        }
        bi.printModel();
    }

    private Resource findVariable(String varName) {
        Query query = QueryFactory.create(
            "PREFIX rdfs: <http://www.w3.org/2000/01/rdf-schema#> " +
            "SELECT ?x " +
            "WHERE {" +
            "    ?x rdfs:label \"" + varName + "\" . " +
            "}");
        QueryExecution qe = QueryExecutionFactory.create(query, model);
        ResultSet results = qe.execSelect();
        QuerySolution qs = results.next();
        Resource x = (Resource) qs.get("x");
        qe.close();
        return x;
    }
    
    private void setBooleanValue(String varName, boolean value) {
        Resource x = findVariable(varName);
        Literal newValue = model.createTypedLiteral(new Boolean(value));
        model.removeAll(x, valueProperty, null);
        model.add(x, valueProperty, newValue);
    }
    
    /**
     * Constructor
     */
    public BayesianInference() {
        model = ModelFactory.createDefaultModel();
        modelReadFile("bayesInf.rdf", model, baseUri);
        valueProperty = model.createProperty(baseUri + "value");
        probProperty = model.createProperty(baseUri + "probability");
    }

    public void printModel() {
        printModel(model);
    }
    
    public float getProbability(String varName) {
        String queryString =
            "PREFIX rdfs: <http://www.w3.org/2000/01/rdf-schema#> " +
            "PREFIX bayes: <" + baseUri + "> " +
            "SELECT ?x ?y " +
            "WHERE {" +
            "    ?x rdfs:label \"" + varName + "\" . " +
            "    ?x bayes:probability ?y . " +
            "}";
        Query query = QueryFactory.create(queryString);
        QueryExecution qe = QueryExecutionFactory.create(query, model);
        ResultSet results = qe.execSelect();
        Literal y = (Literal) results.next().get("y");
        qe.close();
        return y.getFloat();
    }
    
    private final Node_RuleVariable linkVar = new Node_RuleVariable("z", 0);
    private final Node_RuleVariable pbaVar = new Node_RuleVariable("y", 1);
    private final Node_RuleVariable pbnaVar = new Node_RuleVariable("x", 2);
    private final Node_RuleVariable aVar = new Node_RuleVariable("w", 3);
    private final Node_RuleVariable bVar = new Node_RuleVariable("v", 4);
    private final Node_RuleVariable bValueVar = new Node_RuleVariable("u", 5);
    private final Node_RuleVariable probAVar = new Node_RuleVariable("t", 6);

    private final Node probPred = Node.createURI(baseUri + "probability");
    private final Node causePred = Node.createURI(baseUri + "cause");
    private final Node effectPred = Node.createURI(baseUri + "effect");
    private final Node valuePred = Node.createURI(baseUri + "value");
    private final Node p11Pred = Node.createURI(baseUri + "p11");
    private final Node p01Pred = Node.createURI(baseUri + "p01");

    private final ClauseEntry[] updateRuleAntecedents =
        new ClauseEntry[] {
            new TriplePattern(aVar, probPred, probAVar),
            new TriplePattern(linkVar, causePred, aVar),
            new TriplePattern(linkVar, effectPred, bVar),
            new TriplePattern(linkVar, p11Pred, pbaVar),
            new TriplePattern(linkVar, p01Pred, pbnaVar),
            new TriplePattern(bVar, valuePred, bValueVar),
        };
    
    private final Action updateRuleAction =
        new Action() {
            public void run(RuleContext context) {
                Resource aNode = ruleVarToResource(context, model, aVar);
    
                boolean bValue = ruleVarToBoolean(context, model, bValueVar);
                
                float probA = ruleVarToFloat(context, model, probAVar);
                float prBgivenA = ruleVarToFloat(context, model, pbaVar);
                float prBgivenNotA = ruleVarToFloat(context, model, pbnaVar);
    
                float p11 = probA * prBgivenA;
                float p10 = probA * (1.0f - prBgivenA);
                float p01 = (1.0f - probA) * prBgivenNotA;
                float p00 = (1.0f - probA) * (1.0f - prBgivenNotA);
                float probAgivenB = p11 / (p11 + p01);
                float probAgivenNotB = p10 / (p10 + p00);
    
                // update Pr(A) depending on value of B 
                if (bValue) probA = probAgivenB;
                else        probA = probAgivenNotB;
    
                model.removeAll(aNode, probProperty, null);
                model.add(aNode, probProperty,
                    model.createTypedLiteral(new Float(probA)));
            }
        };
        
    public void updateProb() {
        final Rule updateRule =
            buildCustomRule(true,
                    /*
                     * Why doesn't this work?? It should speed this up
                     * because I've pre-selected the variables manually,
                     * but should otherwise work identically.
                    new Node_RuleVariable[] {
                        aVar, bVar, bValueVar, probAVar, linkVar,
                        pbaVar, pbnaVar
                    },
                    */
                    updateRuleAntecedents,
                    updateRuleAction);
        List<Rule> rules = new ArrayList<Rule>();
        rules.add(updateRule);
        Reasoner reasoner = new GenericRuleReasoner(rules);
        reasoner.setDerivationLogging(true);
        model = ModelFactory.createInfModel(reasoner, model);
    }
}