12 files changed, 1710 insertions, 0 deletions
diff --git a/gi/posterior-regularisation/prjava/lib/optimization-0.1.jar b/gi/posterior-regularisation/prjava/lib/optimization-0.1.jar
new file mode 100644
index 00000000..2d14fce7
--- /dev/null
+++ b/gi/posterior-regularisation/prjava/lib/optimization-0.1.jar
diff --git a/gi/posterior-regularisation/prjava/lib/trove-2.0.2.jar b/gi/posterior-regularisation/prjava/lib/trove-2.0.2.jar
new file mode 100644
index 00000000..3e59fbf3
--- /dev/null
+++ b/gi/posterior-regularisation/prjava/lib/trove-2.0.2.jar
diff --git a/gi/posterior-regularisation/prjava/src/data/Corpus.java b/gi/posterior-regularisation/prjava/src/data/Corpus.java
new file mode 100644
index 00000000..f0da0b33
--- /dev/null
+++ b/gi/posterior-regularisation/prjava/src/data/Corpus.java
@@ -0,0 +1,230 @@
+package data;
+
+import java.util.ArrayList;
+import java.util.Arrays;
+import java.util.HashMap;
+import java.util.Scanner;
+
+public class Corpus {
+
+	public static final String alphaFilename="../posdata/corpus.alphabet";
+	public static final String tagalphaFilename="../posdata/corpus.tag.alphabet";
+	
+//	public static final String START_SYM="<s>";
+	public static final String END_SYM="<e>";
+	public static final String NUM_TOK="<NUM>";
+	
+	public static final String UNK_TOK="<unk>";
+	
+	private ArrayList<String[]>sent;
+	private ArrayList<int[]>data;
+	
+	public ArrayList<String[]>tag;
+	public  ArrayList<int[]>tagData;
+	
+	public static boolean convertNumTok=true;
+	
+	private HashMap<String,Integer>freq;
+	public HashMap<String,Integer>vocab;
+	
+	public HashMap<String,Integer>tagVocab;
+	private int tagV;
+	
+	private int V;
+	
+	public static void main(String[] args) {
+		Corpus c=new Corpus("../posdata/en_test.conll");
+		System.out.println(
+			Arrays.toString(c.get(0))	
+		);
+		System.out.println(
+				Arrays.toString(c.getInt(0))	
+			);
+		
+		System.out.println(
+				Arrays.toString(c.get(1))	
+			);
+			System.out.println(
+					Arrays.toString(c.getInt(1))	
+				);
+	}
+
+	public Corpus(String filename,HashMap<String,Integer>dict){
+		V=0;
+		tagV=0;
+		freq=new HashMap<String,Integer>();
+		tagVocab=new HashMap<String,Integer>();
+		vocab=dict;
+		
+		sent=new ArrayList<String[]>();
+		tag=new ArrayList<String[]>();
+		
+		Scanner sc=io.FileUtil.openInFile(filename);
+		ArrayList<String>s=new ArrayList<String>();
+	//	s.add(START_SYM);
+		while(sc.hasNextLine()){
+			String line=sc.nextLine();
+			String toks[]=line.split("\t");
+			if(toks.length<2){
+				s.add(END_SYM);
+				sent.add(s.toArray(new String[0]));
+				s=new ArrayList<String>();
+		//		s.add(START_SYM);
+				continue;
+			}
+			String tok=toks[1].toLowerCase();
+			s.add(tok);
+		}
+		sc.close();
+
+		buildData();
+	}
+	
+	public Corpus(String filename){
+		V=0;
+		freq=new HashMap<String,Integer>();
+		vocab=new HashMap<String,Integer>();
+		tagVocab=new HashMap<String,Integer>();
+		
+		sent=new ArrayList<String[]>();
+		tag=new ArrayList<String[]>();
+		
+		System.out.println("Reading:"+filename);
+		
+		Scanner sc=io.FileUtil.openInFile(filename);
+		ArrayList<String>s=new ArrayList<String>();
+		ArrayList<String>tags=new ArrayList<String>();
+		//s.add(START_SYM);
+		while(sc.hasNextLine()){
+			String line=sc.nextLine();
+			String toks[]=line.split("\t");
+			if(toks.length<2){
+				s.add(END_SYM);
+				tags.add(END_SYM);
+				if(s.size()>2){
+					sent.add(s.toArray(new String[0]));
+					tag.add(tags.toArray(new String [0]));
+				}
+				s=new ArrayList<String>();
+				tags=new ArrayList<String>();
+			//	s.add(START_SYM);
+				continue;
+			}
+			
+			String tok=toks[1].toLowerCase();
+			if(convertNumTok && tok.matches(".*\\d.*")){
+				tok=NUM_TOK;
+			}
+			s.add(tok);
+			
+			if(toks.length>3){
+				tok=toks[3].toLowerCase();
+			}else{
+				tok="_";
+			}
+			tags.add(tok);
+			
+		}
+		sc.close();
+		
+		for(int i=0;i<sent.size();i++){
+			String[]toks=sent.get(i);
+			for(int j=0;j<toks.length;j++){
+				addVocab(toks[j]);
+				addTag(tag.get(i)[j]);
+			}
+		}
+		
+		buildVocab();
+		buildData();
+		System.out.println(data.size()+"sentences, "+vocab.keySet().size()+" word types");
+	}
+
+	public String[] get(int idx){
+		return sent.get(idx);
+	}
+	
+	private void addVocab(String s){
+		Integer integer=freq.get(s);
+		if(integer==null){
+			integer=0;
+		}
+		freq.put(s, integer+1);
+	}
+	
+	public int tokIdx(String tok){
+		Integer integer=vocab.get(tok);
+		if(integer==null){
+			return V;
+		}
+		return integer;
+	}
+	
+	public int tagIdx(String tok){
+		Integer integer=tagVocab.get(tok);
+		if(integer==null){
+			return tagV;
+		}
+		return integer;
+	}
+	
+	private void buildData(){
+		data=new ArrayList<int[]>();
+		for(int i=0;i<sent.size();i++){
+			String s[]=sent.get(i);
+			data.add(new int [s.length]);
+			for(int j=0;j<s.length;j++){
+				data.get(i)[j]=tokIdx(s[j]);
+			}
+		}
+		
+		tagData=new ArrayList<int[]>();
+		for(int i=0;i<tag.size();i++){
+			String s[]=tag.get(i);
+			tagData.add(new int [s.length]);
+			for(int j=0;j<s.length;j++){
+				tagData.get(i)[j]=tagIdx(s[j]);
+			}
+		}
+	}
+	
+	public int [] getInt(int idx){
+		return data.get(idx);
+	}
+	
+	/**
+	 * 
+	 * @return size of vocabulary 
+	 */
+	public int getVocabSize(){
+		return V;
+	}
+	
+	public int [][]getAllData(){
+		return data.toArray(new int [0][]);
+	}
+	
+	public int [][]getTagData(){
+		return tagData.toArray(new int [0][]);
+	}
+	
+	private void buildVocab(){
+		for (String key:freq.keySet()){
+			if(freq.get(key)>2){
+				vocab.put(key, V);
+				V++;
+			}
+		}
+		io.SerializedObjects.writeSerializedObject(vocab, alphaFilename);
+		io.SerializedObjects.writeSerializedObject(tagVocab,tagalphaFilename);
+	}
+
+	private void addTag(String tag){
+		Integer i=tagVocab.get(tag);
+		if(i==null){
+			tagVocab.put(tag, tagV);
+			tagV++;
+		}
+	}
+	
+}
diff --git a/gi/posterior-regularisation/prjava/src/hmm/HMM.java b/gi/posterior-regularisation/prjava/src/hmm/HMM.java
new file mode 100644
index 00000000..1c4d7659
--- /dev/null
+++ b/gi/posterior-regularisation/prjava/src/hmm/HMM.java
@@ -0,0 +1,576 @@
+package hmm;
+
+import java.io.PrintStream;
+import java.util.ArrayList;
+import java.util.Scanner;
+
+public class HMM {
+
+	
+	//trans[i][j]=prob of going FROM i to j
+	double [][]trans;
+	double [][]emit;
+	double []pi;
+	int  [][]data;
+	int [][]tagdata;
+	
+	double logtrans[][];
+	
+	public HMMObjective o;
+	
+	public static void main(String[] args) {
+	
+	}
+	
+	public HMM(int n_state,int n_emit,int [][]data){
+		trans=new double [n_state][n_state];
+		emit=new double[n_state][n_emit];
+		pi=new double [n_state];
+		System.out.println(" random initial parameters");
+		fillRand(trans);
+		fillRand(emit);
+		fillRand(pi);
+
+		this.data=data;
+		
+	}
+	
+	private void fillRand(double [][] a){
+		for(int i=0;i<a.length;i++){
+			for(int j=0;j<a[i].length;j++){
+				a[i][j]=Math.random();
+			}
+			l1normalize(a[i]);
+		}
+	}
+	private void fillRand(double []a){
+		for(int i=0;i<a.length;i++){
+				a[i]=Math.random();
+		}
+		l1normalize(a);
+	}
+	
+	private double loglikely=0;
+	
+	public void EM(){
+		double trans_exp_cnt[][]=new double [trans.length][trans.length];
+		double emit_exp_cnt[][]=new double[trans.length][emit[0].length];
+		double start_exp_cnt[]=new double[trans.length];
+		loglikely=0;
+		
+		//E
+		for(int i=0;i<data.length;i++){
+			
+			double [][][] post=forwardBackward(data[i]);
+			incrementExpCnt(post, data[i], 
+					trans_exp_cnt,
+					emit_exp_cnt,
+					start_exp_cnt);
+			
+			
+			if(i%100==0){
+				System.out.print(".");
+			}
+			if(i%1000==0){
+				System.out.println(i);
+			}
+			
+		}
+		System.out.println("Log likelihood: "+loglikely);
+		
+		//M
+		addOneSmooth(emit_exp_cnt);
+		for(int i=0;i<trans.length;i++){
+		
+			//transition probs
+			double sum=0;
+			for(int j=0;j<trans.length;j++){
+				sum+=trans_exp_cnt[i][j];
+			}
+			//avoid NAN
+			if(sum==0){
+				sum=1;
+			}
+			for(int j=0;j<trans[i].length;j++){
+				trans[i][j]=trans_exp_cnt[i][j]/sum;
+			}
+			
+			//emission probs
+			
+			sum=0;
+			for(int j=0;j<emit[i].length;j++){
+				sum+=emit_exp_cnt[i][j];
+			}
+			//avoid NAN
+			if(sum==0){
+				sum=1;
+			}
+			for(int j=0;j<emit[i].length;j++){
+				emit[i][j]=emit_exp_cnt[i][j]/sum;
+			}
+			
+			
+			//initial probs
+			for(int j=0;j<pi.length;j++){
+				pi[j]=start_exp_cnt[j];
+			}
+			l1normalize(pi);
+		}
+	}
+	
+	private double [][][]forwardBackward(int [] seq){
+		double a[][]=new double [seq.length][trans.length];
+		double b[][]=new double [seq.length][trans.length];
+		
+		int len=seq.length;
+		//initialize the first step
+		for(int i=0;i<trans.length;i++){
+			a[0][i]=emit[i][seq[0]]*pi[i];
+			b[len-1][i]=1;
+		}
+		
+		//log of denominator for likelyhood
+		double c=Math.log(l1norm(a[0]));
+		
+		l1normalize(a[0]);
+		l1normalize(b[len-1]);
+		
+		
+		
+		//forward
+		for(int n=1;n<len;n++){
+			for(int i=0;i<trans.length;i++){
+				for(int j=0;j<trans.length;j++){
+					a[n][i]+=trans[j][i]*a[n-1][j];
+				}
+				a[n][i]*=emit[i][seq[n]];
+			}
+			c+=Math.log(l1norm(a[n]));
+			l1normalize(a[n]);
+		}
+		
+		loglikely+=c;
+		
+		//backward
+		for(int n=len-2;n>=0;n--){
+			for(int i=0;i<trans.length;i++){
+				for(int j=0;j<trans.length;j++){
+					b[n][i]+=trans[i][j]*b[n+1][j]*emit[j][seq[n+1]];
+				}
+			}
+			l1normalize(b[n]);
+		}
+		
+		
+		//expected transition 
+		double p[][][]=new double [seq.length][trans.length][trans.length];
+		for(int n=0;n<len-1;n++){
+			for(int i=0;i<trans.length;i++){
+				for(int j=0;j<trans.length;j++){
+					p[n][i][j]=a[n][i]*trans[i][j]*emit[j][seq[n+1]]*b[n+1][j];
+					
+				}
+			}
+
+			l1normalize(p[n]);
+		}
+		return p;
+	}
+	
+	private void incrementExpCnt(
+			double post[][][],int [] seq, 
+			double trans_exp_cnt[][],
+			double emit_exp_cnt[][],
+			double start_exp_cnt[])
+	{
+		
+		for(int n=0;n<post.length;n++){
+			for(int i=0;i<trans.length;i++){
+				double py=0;
+				for(int j=0;j<trans.length;j++){
+					py+=post[n][i][j];
+					trans_exp_cnt[i][j]+=post[n][i][j];
+				}
+
+				emit_exp_cnt[i][seq[n]]+=py;				
+				
+			}
+		}
+		
+		//the first state
+		for(int i=0;i<trans.length;i++){
+			double py=0;
+			for(int j=0;j<trans.length;j++){
+				py+=post[0][i][j];
+			}
+			start_exp_cnt[i]+=py;	
+		}
+		
+		
+		//the last state
+		int len=post.length;
+		for(int i=0;i<trans.length;i++){
+			double py=0;
+			for(int j=0;j<trans.length;j++){
+				py+=post[len-2][j][i];
+			}
+			emit_exp_cnt[i][seq[len-1]]+=py;	
+		}
+	}
+	
+	public void l1normalize(double [] a){
+		double sum=0;
+		for(int i=0;i<a.length;i++){
+			sum+=a[i];
+		}
+		if(sum==0){
+			return ;
+		}
+		for(int i=0;i<a.length;i++){
+			a[i]/=sum;
+		}
+	}
+	
+	public  void l1normalize(double [][] a){
+		double sum=0;
+		for(int i=0;i<a.length;i++){
+			for(int j=0;j<a[i].length;j++){
+				sum+=a[i][j];
+			}
+		}
+		if(sum==0){
+			return;
+		}
+		for(int i=0;i<a.length;i++){
+			for(int j=0;j<a[i].length;j++){
+				a[i][j]/=sum;
+			}
+		}
+	}
+	
+	public void writeModel(String modelFilename){
+		PrintStream ps=io.FileUtil.openOutFile(modelFilename);
+		ps.println(trans.length);
+		ps.println("Initial Probabilities:");
+		for(int i=0;i<pi.length;i++){
+			ps.print(pi[i]+"\t");
+		}
+		ps.println();
+		ps.println("Transition Probabilities:");
+		for(int i=0;i<trans.length;i++){
+			for(int j=0;j<trans[i].length;j++){
+				ps.print(trans[i][j]+"\t");
+			}
+			ps.println();
+		}
+		ps.println("Emission Probabilities:");
+		ps.println(emit[0].length);
+		for(int i=0;i<trans.length;i++){
+			for(int j=0;j<emit[i].length;j++){
+				ps.println(emit[i][j]);
+			}
+			ps.println();
+		}
+		ps.close();
+	}
+	
+	public HMM(){
+	
+	}
+	
+	public void readModel(String modelFilename){
+		Scanner sc=io.FileUtil.openInFile(modelFilename);
+	
+		int n_state=sc.nextInt();
+		sc.nextLine();
+		sc.nextLine();
+		pi=new double [n_state];
+		for(int i=0;i<n_state;i++){
+			pi[i]=sc.nextDouble();
+		}
+		sc.nextLine();
+		sc.nextLine();
+		trans=new double[n_state][n_state];
+		for(int i=0;i<trans.length;i++){
+			for(int j=0;j<trans[i].length;j++){
+				trans[i][j]=sc.nextDouble();
+			}
+		}
+		sc.nextLine();
+		sc.nextLine();
+		
+		int n_obs=sc.nextInt();
+		emit=new double[n_state][n_obs];
+		for(int i=0;i<trans.length;i++){
+			for(int j=0;j<emit[i].length;j++){
+				emit[i][j]=sc.nextDouble();
+			}
+		}
+		sc.close();
+	}
+	
+	public int []viterbi(int [] seq){
+		double [][]p=new double [seq.length][trans.length];
+		int backp[][]=new int [seq.length][trans.length];
+		
+		for(int i=0;i<trans.length;i++){
+			p[0][i]=Math.log(emit[i][seq[0]]*pi[i]);
+		}
+		
+		double a[][]=logtrans;
+		if(logtrans==null){
+			a=new double [trans.length][trans.length];
+			for(int i=0;i<trans.length;i++){
+				for(int j=0;j<trans.length;j++){
+					a[i][j]=Math.log(trans[i][j]);
+				}
+			}
+			logtrans=a;
+		}
+		
+		double maxprob=0;
+		for(int n=1;n<seq.length;n++){
+			for(int i=0;i<trans.length;i++){
+				maxprob=p[n-1][0]+a[0][i];
+				backp[n][i]=0;
+				for(int j=1;j<trans.length;j++){
+					double prob=p[n-1][j]+a[j][i];
+					if(maxprob<prob){
+						backp[n][i]=j;
+						maxprob=prob;
+					}
+				}
+				p[n][i]=maxprob+Math.log(emit[i][seq[n]]);
+			}
+		}
+		
+		maxprob=p[seq.length-1][0];
+		int maxIdx=0;
+		for(int i=1;i<trans.length;i++){
+			if(p[seq.length-1][i]>maxprob){
+				maxprob=p[seq.length-1][i];
+				maxIdx=i;
+			}
+		}
+		int ans[]=new int [seq.length];
+		ans[seq.length-1]=maxIdx;
+		for(int i=seq.length-2;i>=0;i--){
+			ans[i]=backp[i+1][ans[i+1]];
+		}
+		return ans;
+	}
+	
+	public double l1norm(double a[]){
+		double norm=0;
+		for(int i=0;i<a.length;i++){
+			norm += a[i];
+		}
+		return norm;
+	}
+	
+	public double [][]getEmitProb(){
+		return emit;
+	}
+	
+	public int [] sample(int terminalSym){
+		ArrayList<Integer > s=new ArrayList<Integer>();
+		int state=sample(pi);
+		int sym=sample(emit[state]);
+		while(sym!=terminalSym){
+			s.add(sym);
+			state=sample(trans[state]);
+			sym=sample(emit[state]);
+		}
+		
+		int ans[]=new int [s.size()];
+		for(int i=0;i<ans.length;i++){
+			ans[i]=s.get(i);
+		}
+		return ans;
+	}
+	
+	public int sample(double p[]){
+		double r=Math.random();
+		double sum=0;
+		for(int i=0;i<p.length;i++){
+			sum+=p[i];
+			if(sum>=r){
+				return i;
+			}
+		}
+		return p.length-1;
+	}
+	
+	public void train(int tagdata[][]){
+		double trans_exp_cnt[][]=new double [trans.length][trans.length];
+		double emit_exp_cnt[][]=new double[trans.length][emit[0].length];
+		double start_exp_cnt[]=new double[trans.length];
+		
+		for(int i=0;i<tagdata.length;i++){
+			start_exp_cnt[tagdata[i][0]]++;
+			
+			for(int j=0;j<tagdata[i].length;j++){
+				if(j+1<tagdata[i].length){
+					trans_exp_cnt[ tagdata[i][j] ] [ tagdata[i][j+1] ]++;
+				}
+				emit_exp_cnt[tagdata[i][j]][data[i][j]]++;
+			}
+			
+		}
+		
+		//M
+		addOneSmooth(emit_exp_cnt);
+		for(int i=0;i<trans.length;i++){
+		
+			//transition probs
+			double sum=0;
+			for(int j=0;j<trans.length;j++){
+				sum+=trans_exp_cnt[i][j];
+			}
+			if(sum==0){
+				sum=1;
+			}
+			for(int j=0;j<trans[i].length;j++){
+				trans[i][j]=trans_exp_cnt[i][j]/sum;
+			}
+			
+			//emission probs
+
+			sum=0;
+			for(int j=0;j<emit[i].length;j++){
+				sum+=emit_exp_cnt[i][j];
+			}
+			if(sum==0){
+				sum=1;
+			}
+			for(int j=0;j<emit[i].length;j++){
+				emit[i][j]=emit_exp_cnt[i][j]/sum;
+			}
+
+			
+			//initial probs
+			for(int j=0;j<pi.length;j++){
+				pi[j]=start_exp_cnt[j];
+			}
+			l1normalize(pi);
+		}
+	}
+	
+	private void addOneSmooth(double a[][]){
+		for(int i=0;i<a.length;i++){
+			for(int j=0;j<a[i].length;j++){
+				a[i][j]+=0.01;
+			}
+			//l1normalize(a[i]);
+		}
+	}
+	
+	public void PREM(){
+		
+		o.optimizeWithProjectedGradientDescent();
+		
+		double trans_exp_cnt[][]=new double [trans.length][trans.length];
+		double emit_exp_cnt[][]=new double[trans.length][emit[0].length];
+		double start_exp_cnt[]=new double[trans.length];
+		
+		o.loglikelihood=0;
+		//E
+		for(int sentNum=0;sentNum<data.length;sentNum++){
+			
+			double [][][] post=o.forwardBackward(sentNum);
+			incrementExpCnt(post, data[sentNum], 
+					trans_exp_cnt,
+					emit_exp_cnt,
+					start_exp_cnt);
+			
+			
+			if(sentNum%100==0){
+				System.out.print(".");
+			}
+			if(sentNum%1000==0){
+				System.out.println(sentNum);
+			}
+			
+		}
+		
+		System.out.println("Log likelihood: "+o.getValue());
+		
+		//M
+		addOneSmooth(emit_exp_cnt);
+		for(int i=0;i<trans.length;i++){
+		
+			//transition probs
+			double sum=0;
+			for(int j=0;j<trans.length;j++){
+				sum+=trans_exp_cnt[i][j];
+			}
+			//avoid NAN
+			if(sum==0){
+				sum=1;
+			}
+			for(int j=0;j<trans[i].length;j++){
+				trans[i][j]=trans_exp_cnt[i][j]/sum;
+			}
+			
+			//emission probs
+			
+			sum=0;
+			for(int j=0;j<emit[i].length;j++){
+				sum+=emit_exp_cnt[i][j];
+			}
+			//avoid NAN
+			if(sum==0){
+				sum=1;
+			}
+			for(int j=0;j<emit[i].length;j++){
+				emit[i][j]=emit_exp_cnt[i][j]/sum;
+			}
+			
+			
+			//initial probs
+			for(int j=0;j<pi.length;j++){
+				pi[j]=start_exp_cnt[j];
+			}
+			l1normalize(pi);
+		}
+		
+	}
+	
+	public void computeMaxwt(double[][]maxwt, int[][] d){
+
+		for(int sentNum=0;sentNum<d.length;sentNum++){
+			double post[][][]=forwardBackward(d[sentNum]);
+			
+			for(int n=0;n<post.length;n++){
+				for(int i=0;i<trans.length;i++){
+					double py=0;
+					for(int j=0;j<trans.length;j++){
+						py+=post[n][i][j];
+					}
+
+					if(py>maxwt[i][d[sentNum][n]]){
+						maxwt[i][d[sentNum][n]]=py;
+					}
+					
+				}
+			}
+			
+			//the last state
+			int len=post.length;
+			for(int i=0;i<trans.length;i++){
+				double py=0;
+				for(int j=0;j<trans.length;j++){
+					py+=post[len-2][j][i];
+				}
+				
+				if(py>maxwt[i][d[sentNum][len-1]]){
+					maxwt[i][d[sentNum][len-1]]=py;
+				}
+				
+			}
+			
+		}
+	
+	}
+	
+}//end of class
diff --git a/gi/posterior-regularisation/prjava/src/hmm/HMMObjective.java b/gi/posterior-regularisation/prjava/src/hmm/HMMObjective.java
new file mode 100644
index 00000000..551210c0
--- /dev/null
+++ b/gi/posterior-regularisation/prjava/src/hmm/HMMObjective.java
@@ -0,0 +1,348 @@
+package hmm;
+
+import gnu.trove.TIntArrayList;
+import optimization.gradientBasedMethods.ProjectedGradientDescent;
+import optimization.gradientBasedMethods.ProjectedObjective;
+import optimization.gradientBasedMethods.stats.OptimizerStats;
+import optimization.linesearch.ArmijoLineSearchMinimizationAlongProjectionArc;
+import optimization.linesearch.InterpolationPickFirstStep;
+import optimization.linesearch.LineSearchMethod;
+import optimization.projections.SimplexProjection;
+import optimization.stopCriteria.CompositeStopingCriteria;
+import optimization.stopCriteria.ProjectedGradientL2Norm;
+import optimization.stopCriteria.StopingCriteria;
+import optimization.stopCriteria.ValueDifference;
+
+public class HMMObjective extends ProjectedObjective{
+
+	
+	private static final double GRAD_DIFF = 3;
+	public static double INIT_STEP_SIZE=10;
+	public static double VAL_DIFF=2000;
+	
+	private HMM hmm;
+	double[] newPoint  ;
+	
+	//posterior[sent num][tok num][tag]=index into lambda
+	private int posteriorMap[][][];
+	//projection[word][tag].get(occurence)=index into lambda
+	private TIntArrayList projectionMap[][];
+
+	//Size of the simplex
+	public double scale=10;
+	private SimplexProjection projection;
+	
+	private int wordFreq[];
+	private static int MIN_FREQ=3;
+	private int numWordsToProject=0;
+	
+	private int n_param;
+	
+	public  double loglikelihood;
+	
+	public HMMObjective(HMM h){
+		hmm=h;
+		
+		countWords();
+		buildMap();
+
+		gradient=new double [n_param];
+		projection = new SimplexProjection(scale);
+		newPoint  = new double[n_param];
+		setInitialParameters(new double[n_param]);
+		
+	}
+	
+	/**@brief counts word frequency in the corpus
+	 * 
+	 */
+	private void countWords(){
+		wordFreq=new int [hmm.emit[0].length];
+		for(int i=0;i<hmm.data.length;i++){
+			for(int j=0;j<hmm.data[i].length;j++){
+				wordFreq[hmm.data[i][j]]++;
+			}
+		}
+	}
+	
+	/**@brief build posterior and projection indices
+	 * 
+	 */
+	private void buildMap(){
+		//number of sentences hidden states and words
+		int n_states=hmm.trans.length;
+		int n_words=hmm.emit[0].length;
+		int n_sents=hmm.data.length;
+		
+		n_param=0;
+		posteriorMap=new int[n_sents][][];
+		projectionMap=new TIntArrayList[n_words][];
+		for(int sentNum=0;sentNum<n_sents;sentNum++){
+			int [] data=hmm.data[sentNum];
+			posteriorMap[sentNum]=new int[data.length][n_states];
+			numWordsToProject=0;
+			for(int i=0;i<data.length;i++){
+				int word=data[i];
+				for(int state=0;state<n_states;state++){
+					if(wordFreq[word]>MIN_FREQ){
+						if(projectionMap[word]==null){
+							projectionMap[word]=new TIntArrayList[n_states];
+						}
+						
+						posteriorMap[sentNum][i][state]=n_param;
+						if(projectionMap[word][state]==null){
+							projectionMap[word][state]=new TIntArrayList();
+							numWordsToProject++;
+						}
+						projectionMap[word][state].add(n_param);
+						n_param++;
+					}else{
+						
+						posteriorMap[sentNum][i][state]=-1;
+					}
+				}
+			}
+		}
+	}
+	
+	@Override
+	public double[] projectPoint(double[] point) {
+		// TODO Auto-generated method stub
+		for(int i=0;i<projectionMap.length;i++){
+			
+			if(projectionMap[i]==null){
+				//this word is not constrained
+				continue;
+			}
+			
+			for(int j=0;j<projectionMap[i].length;j++){
+				TIntArrayList instances=projectionMap[i][j];
+				double[] toProject = new double[instances.size()];
+				
+				for (int k = 0; k < toProject.length; k++) {
+					//	System.out.print(instances.get(k) + " ");
+						toProject[k] = point[instances.get(k)];
+				}
+				
+				projection.project(toProject);
+				for (int k = 0; k < toProject.length; k++) {
+					newPoint[instances.get(k)]=toProject[k];
+				}
+			}
+		}
+		return newPoint;
+	}
+
+	@Override
+	public double[] getGradient() {
+		// TODO Auto-generated method stub
+		gradientCalls++;
+		return gradient;
+	}
+
+	@Override
+	public double getValue() {
+		// TODO Auto-generated method stub
+		functionCalls++;
+		return loglikelihood;
+	}
+	
+
+	@Override
+	public String toString() {
+		// TODO Auto-generated method stub
+		StringBuffer sb = new StringBuffer();
+		for (int i = 0; i < parameters.length; i++) {
+			sb.append(parameters[i]+" ");
+			if(i%100==0){
+				sb.append("\n");
+			}
+		}
+		sb.append("\n");
+		/*
+		for (int i = 0; i < gradient.length; i++) {
+			sb.append(gradient[i]+" ");
+			if(i%100==0){
+				sb.append("\n");
+			}
+		}
+		sb.append("\n");
+		*/
+		return sb.toString();
+	}
+
+	
+	/**
+	 * @param seq
+	 * @return posterior probability of each transition
+	 */
+	public double [][][]forwardBackward(int sentNum){
+		int [] seq=hmm.data[sentNum];
+		int n_states=hmm.trans.length;
+		double a[][]=new double [seq.length][n_states];
+		double b[][]=new double [seq.length][n_states];
+		
+		int len=seq.length;
+		
+		boolean  constrained=
+			(projectionMap[seq[0]]!=null);
+
+		//initialize the first step
+		for(int i=0;i<n_states;i++){
+			a[0][i]=hmm.emit[i][seq[0]]*hmm.pi[i];
+			if(constrained){
+				a[0][i]*=
+					Math.exp(- parameters[ posteriorMap[sentNum][0][i] ] );
+			}
+			b[len-1][i]=1;
+		}
+		
+		loglikelihood+=Math.log(hmm.l1norm(a[0]));		
+		hmm.l1normalize(a[0]);
+		hmm.l1normalize(b[len-1]);
+		
+		//forward
+		for(int n=1;n<len;n++){
+			
+			constrained=
+				(projectionMap[seq[n]]!=null);
+			
+			for(int i=0;i<n_states;i++){
+				for(int j=0;j<n_states;j++){
+					a[n][i]+=hmm.trans[j][i]*a[n-1][j];
+				}
+				a[n][i]*=hmm.emit[i][seq[n]];
+				
+				if(constrained){
+					a[n][i]*=
+						Math.exp(- parameters[ posteriorMap[sentNum][n][i] ] );
+				}
+				
+			}
+			loglikelihood+=Math.log(hmm.l1norm(a[n]));
+			hmm.l1normalize(a[n]);
+		}
+		
+		//temp variable for e^{-\lambda}
+		double factor=1;
+		//backward
+		for(int n=len-2;n>=0;n--){
+			
+			constrained=
+				(projectionMap[seq[n+1]]!=null);
+			
+			for(int i=0;i<n_states;i++){
+				for(int j=0;j<n_states;j++){
+					
+					if(constrained){
+						factor=
+							Math.exp(- parameters[ posteriorMap[sentNum][n+1][j] ] );
+					}else{
+						factor=1;
+					}
+					
+					b[n][i]+=hmm.trans[i][j]*b[n+1][j]*hmm.emit[j][seq[n+1]]*factor;
+					
+				}
+			}
+			hmm.l1normalize(b[n]);
+		}
+		
+		//expected transition 
+		double p[][][]=new double [seq.length][n_states][n_states];
+		for(int n=0;n<len-1;n++){
+			
+			constrained=
+				(projectionMap[seq[n+1]]!=null);
+			
+			for(int i=0;i<n_states;i++){
+				for(int j=0;j<n_states;j++){
+					
+					if(constrained){
+						factor=
+							Math.exp(- parameters[ posteriorMap[sentNum][n+1][j] ] );
+					}else{
+						factor=1;
+					}
+					
+					p[n][i][j]=a[n][i]*hmm.trans[i][j]*
+						hmm.emit[j][seq[n+1]]*b[n+1][j]*factor;
+					
+				}
+			}
+
+			hmm.l1normalize(p[n]);
+		}
+		return p;
+	}
+
+	public void optimizeWithProjectedGradientDescent(){
+		LineSearchMethod ls =
+			new ArmijoLineSearchMinimizationAlongProjectionArc
+				(new InterpolationPickFirstStep(INIT_STEP_SIZE));
+		
+		OptimizerStats stats = new OptimizerStats();
+		
+		
+		ProjectedGradientDescent optimizer = new ProjectedGradientDescent(ls);
+		StopingCriteria stopGrad = new ProjectedGradientL2Norm(GRAD_DIFF);
+		StopingCriteria stopValue = new ValueDifference(VAL_DIFF);
+		CompositeStopingCriteria compositeStop = new CompositeStopingCriteria();
+		compositeStop.add(stopGrad);
+		compositeStop.add(stopValue);
+		
+		optimizer.setMaxIterations(10);
+		updateFunction();
+		boolean succed = optimizer.optimize(this,stats,compositeStop);
+		System.out.println("Ended optimzation Projected Gradient Descent\n" + stats.prettyPrint(1));
+		if(succed){
+			System.out.println("Ended optimization in " + optimizer.getCurrentIteration());
+		}else{
+			System.out.println("Failed to optimize");
+		}
+	}
+	
+	@Override
+	public void setParameters(double[] params) {
+		super.setParameters(params);
+		updateFunction();
+	}
+	
+	private void updateFunction(){
+		
+		updateCalls++;
+		loglikelihood=0;
+	
+		for(int sentNum=0;sentNum<hmm.data.length;sentNum++){
+			double [][][]p=forwardBackward(sentNum);
+			
+			for(int n=0;n<p.length-1;n++){
+				for(int i=0;i<p[n].length;i++){
+					if(projectionMap[hmm.data[sentNum][n]]!=null){
+						double posterior=0;
+						for(int j=0;j<p[n][i].length;j++){
+							posterior+=p[n][i][j];
+						}
+						gradient[posteriorMap[sentNum][n][i]]=-posterior;
+					}
+				}
+			}
+			
+			//the last state
+			int n=p.length-2;
+			for(int i=0;i<p[n].length;i++){
+				if(projectionMap[hmm.data[sentNum][n+1]]!=null){
+					
+					double posterior=0;
+					for(int j=0;j<p[n].length;j++){
+						posterior+=p[n][j][i];
+					}
+					gradient[posteriorMap[sentNum][n+1][i]]=-posterior;
+				
+				}
+			}	
+		}
+		
+	}
+	
+}
diff --git a/gi/posterior-regularisation/prjava/src/hmm/POS.java b/gi/posterior-regularisation/prjava/src/hmm/POS.java
new file mode 100644
index 00000000..722d38e2
--- /dev/null
+++ b/gi/posterior-regularisation/prjava/src/hmm/POS.java
@@ -0,0 +1,126 @@
+package hmm;
+
+import java.io.PrintStream;
+import java.util.HashMap;
+
+import data.Corpus;
+
+public class POS {
+
+	//public String trainFilename="../posdata/en_train.conll";
+	//public static String trainFilename="../posdata/small_train.txt";
+	public static String trainFilename="../posdata/en_test.conll";
+//	public static String trainFilename="../posdata/trial1.txt";
+	
+	public static String testFilename="../posdata/en_test.conll";
+	//public static String testFilename="../posdata/trial1.txt";
+	
+	public static String predFilename="../posdata/en_test.predict.conll";
+	public static String modelFilename="../posdata/posModel.out";
+	public static final int ITER=20;
+	public static final int N_STATE=30;
+	
+	public static void main(String[] args) {
+		//POS p=new POS();
+		//POS p=new POS(true);
+		PRPOS();
+	}
+
+	
+	public POS(){
+		Corpus c= new Corpus(trainFilename);
+		//size of vocabulary +1 for unknown tokens
+		HMM hmm =new HMM(N_STATE, c.getVocabSize()+1,c.getAllData());
+		for(int i=0;i<ITER;i++){
+			System.out.println("Iter"+i);
+			hmm.EM();
+			if((i+1)%10==0){
+				hmm.writeModel(modelFilename+i);
+			}
+		}
+
+		hmm.writeModel(modelFilename);
+		
+		Corpus test=new Corpus(testFilename,c.vocab);
+		
+		PrintStream ps= io.FileUtil.openOutFile(predFilename);
+		
+		int [][]data=test.getAllData();
+		for(int i=0;i<data.length;i++){
+			int []tag=hmm.viterbi(data[i]);
+			String sent[]=test.get(i);
+			for(int j=0;j<data[i].length;j++){
+				ps.println(sent[j]+"\t"+tag[j]);
+			}
+			ps.println();
+		}
+		ps.close();
+	}
+	
+	//POS induction with L1/Linf constraints
+	public static void PRPOS(){
+		Corpus c= new Corpus(trainFilename);
+		//size of vocabulary +1 for unknown tokens
+		HMM hmm =new HMM(N_STATE, c.getVocabSize()+1,c.getAllData());
+		hmm.o=new HMMObjective(hmm);
+		for(int i=0;i<ITER;i++){
+			System.out.println("Iter: "+i);
+			hmm.PREM();
+			if((i+1)%10==0){
+				hmm.writeModel(modelFilename+i);
+			}
+		}
+
+		hmm.writeModel(modelFilename);
+		
+		Corpus test=new Corpus(testFilename,c.vocab);
+		
+		PrintStream ps= io.FileUtil.openOutFile(predFilename);
+		
+		int [][]data=test.getAllData();
+		for(int i=0;i<data.length;i++){
+			int []tag=hmm.viterbi(data[i]);
+			String sent[]=test.get(i);
+			for(int j=0;j<data[i].length;j++){
+				ps.println(sent[j]+"\t"+tag[j]);
+			}
+			ps.println();
+		}
+		ps.close();
+	}
+	
+	
+	public POS(boolean supervised){
+		Corpus c= new Corpus(trainFilename);
+		//size of vocabulary +1 for unknown tokens
+		HMM hmm =new HMM(c.tagVocab.size() , c.getVocabSize()+1,c.getAllData());
+		hmm.train(c.getTagData());
+
+		hmm.writeModel(modelFilename);
+		
+		Corpus test=new Corpus(testFilename,c.vocab);
+		
+		HashMap<String, Integer>tagVocab=
+			(HashMap<String, Integer>) io.SerializedObjects.readSerializedObject(Corpus.tagalphaFilename);
+		String [] tagdict=new String [tagVocab.size()+1];
+		for(String key:tagVocab.keySet()){
+			tagdict[tagVocab.get(key)]=key;
+		}
+		tagdict[tagdict.length-1]=Corpus.UNK_TOK;
+		
+		System.out.println(c.vocab.get("<e>"));
+		
+		PrintStream ps= io.FileUtil.openOutFile(predFilename);
+		
+		int [][]data=test.getAllData();
+		for(int i=0;i<data.length;i++){
+			int []tag=hmm.viterbi(data[i]);
+			String sent[]=test.get(i);
+			for(int j=0;j<data[i].length;j++){
+				ps.println(sent[j]+"\t"+tagdict[tag[j]]);
+			}
+			ps.println();
+		}
+		ps.close();
+	}
+}
diff --git a/gi/posterior-regularisation/prjava/src/io/FileUtil.java b/gi/posterior-regularisation/prjava/src/io/FileUtil.java
new file mode 100644
index 00000000..bfa58213
--- /dev/null
+++ b/gi/posterior-regularisation/prjava/src/io/FileUtil.java
@@ -0,0 +1,37 @@
+package io;
+import java.util.*;
+import java.io.*;
+public class FileUtil {
+	public static Scanner openInFile(String filename){
+		Scanner localsc=null;
+		try
+		{
+			localsc=new Scanner (new FileInputStream(filename));
+
+		}catch(IOException ioe){
+			System.out.println(ioe.getMessage());
+		}
+		return localsc;
+	}
+	public static PrintStream  openOutFile(String filename){
+		PrintStream localps=null;
+		try
+		{
+			localps=new PrintStream (new FileOutputStream(filename));
+
+		}catch(IOException ioe){
+			System.out.println(ioe.getMessage());
+		}
+		return localps;
+	}
+	public static FileInputStream openInputStream(String infilename){
+		FileInputStream fis=null;
+		try {
+			fis =(new FileInputStream(infilename));
+			
+		} catch (IOException ioe) {
+			System.out.println(ioe.getMessage());
+		}
+		return fis;
+	}
+}
diff --git a/gi/posterior-regularisation/prjava/src/io/SerializedObjects.java b/gi/posterior-regularisation/prjava/src/io/SerializedObjects.java
new file mode 100644
index 00000000..d1631b51
--- /dev/null
+++ b/gi/posterior-regularisation/prjava/src/io/SerializedObjects.java
@@ -0,0 +1,83 @@
+package io;
+
+
+
+import java.io.BufferedInputStream;
+import java.io.BufferedOutputStream;
+import java.io.FileInputStream;
+import java.io.FileOutputStream;
+import java.io.IOException;
+import java.io.InputStream;
+import java.io.ObjectInput;
+import java.io.ObjectInputStream;
+import java.io.ObjectOutput;
+import java.io.ObjectOutputStream;
+import java.io.OutputStream;
+
+public class SerializedObjects
+{
+	public static void writeSerializedObject(Object object, String outFile)
+	{
+		ObjectOutput output = null;
+	    try{
+	      //use buffering
+	      OutputStream file = new FileOutputStream(outFile);
+	      OutputStream buffer = new BufferedOutputStream( file );
+	      output = new ObjectOutputStream( buffer );
+	      output.writeObject(object);
+	      buffer.close();
+	      file.close();
+	    }
+	    catch(IOException ex){
+	    	ex.printStackTrace();
+	    }
+	    finally{
+	      try {
+	        if (output != null) {
+	          //flush and close "output" and its underlying streams
+	          output.close();
+	        }
+	      }
+	      catch (IOException ex ){
+	    	  ex.printStackTrace();
+	      }
+	    }
+	}
+	
+	public static Object readSerializedObject(String inputFile)
+	{
+		ObjectInput input = null;
+	    Object recoveredObject=null;
+		try{
+	      //use buffering
+	      InputStream file = new FileInputStream(inputFile);
+	      InputStream buffer = new BufferedInputStream(file);
+	      input = new ObjectInputStream(buffer);
+	      //deserialize the List
+	      recoveredObject = input.readObject();
+	    }
+	    catch(IOException ex){
+	    	ex.printStackTrace();
+	    }
+	    catch (ClassNotFoundException ex){
+	    	ex.printStackTrace();
+	    }
+	    catch(Exception ex)
+	    {
+	    	ex.printStackTrace();
+	    }
+	    finally{
+	      try {
+	        if ( input != null ) {
+	          //close "input" and its underlying streams
+	          input.close();
+	        }
+	      }
+	      catch (IOException ex){
+	    	  ex.printStackTrace();
+	      }
+	    }
+	    return recoveredObject;
+	 }
+	
+}
+\ No newline at end of file
diff --git a/gi/posterior-regularisation/prjava/src/test/CorpusTest.java b/gi/posterior-regularisation/prjava/src/test/CorpusTest.java
new file mode 100644
index 00000000..b4c3041f
--- /dev/null
+++ b/gi/posterior-regularisation/prjava/src/test/CorpusTest.java
@@ -0,0 +1,60 @@
+package test;
+
+import java.util.Arrays;
+import java.util.HashMap;
+
+import data.Corpus;
+import hmm.POS;
+
+public class CorpusTest {
+
+	public static void main(String[] args) {
+		Corpus c=new Corpus(POS.trainFilename);
+
+		
+		int idx=30;
+		
+		
+		HashMap<String, Integer>vocab=
+			(HashMap<String, Integer>) io.SerializedObjects.readSerializedObject(Corpus.alphaFilename);
+		
+		HashMap<String, Integer>tagVocab=
+			(HashMap<String, Integer>) io.SerializedObjects.readSerializedObject(Corpus.tagalphaFilename);
+		
+		
+		String [] dict=new String [vocab.size()+1];
+		for(String key:vocab.keySet()){
+			dict[vocab.get(key)]=key;
+		}
+		dict[dict.length-1]=Corpus.UNK_TOK;
+		
+		String [] tagdict=new String [tagVocab.size()+1];
+		for(String key:tagVocab.keySet()){
+			tagdict[tagVocab.get(key)]=key;
+		}
+		tagdict[tagdict.length-1]=Corpus.UNK_TOK;
+		
+		String[] sent=c.get(idx);
+		int []data=c.getInt(idx);
+		
+		
+		String []roundtrip=new String [sent.length];
+		for(int i=0;i<sent.length;i++){
+			roundtrip[i]=dict[data[i]];
+		}
+		System.out.println(Arrays.toString(sent));
+		System.out.println(Arrays.toString(roundtrip));
+		
+		sent=c.tag.get(idx);
+		data=c.tagData.get(idx);
+		
+		
+		roundtrip=new String [sent.length];
+		for(int i=0;i<sent.length;i++){
+			roundtrip[i]=tagdict[data[i]];
+		}
+		System.out.println(Arrays.toString(sent));
+		System.out.println(Arrays.toString(roundtrip));
+	}
+
+}
diff --git a/gi/posterior-regularisation/prjava/src/test/HMMModelStats.java b/gi/posterior-regularisation/prjava/src/test/HMMModelStats.java
new file mode 100644
index 00000000..26d7abec
--- /dev/null
+++ b/gi/posterior-regularisation/prjava/src/test/HMMModelStats.java
@@ -0,0 +1,96 @@
+package test;
+
+import hmm.HMM;
+import hmm.POS;
+
+import java.io.PrintStream;
+import java.util.ArrayList;
+import java.util.Collections;
+import java.util.HashMap;
+
+import data.Corpus;
+
+public class HMMModelStats {
+
+	public static String modelFilename="../posdata/posModel.out";
+	public static String alphaFilename="../posdata/corpus.alphabet";
+	public static String statsFilename="../posdata/model.stats";
+
+	public static final int NUM_WORD=50;
+	
+	public static String testFilename="../posdata/en_test.conll";
+	
+	public static double [][]maxwt;
+	
+	public static void main(String[] args) {
+		HashMap<String, Integer>vocab=
+			(HashMap<String, Integer>) io.SerializedObjects.readSerializedObject(alphaFilename);
+		
+		Corpus test=new Corpus(testFilename,vocab);
+		
+		String [] dict=new String [vocab.size()+1];
+		for(String key:vocab.keySet()){
+			dict[vocab.get(key)]=key;
+		}
+		dict[dict.length-1]=Corpus.UNK_TOK;
+		
+		HMM hmm=new HMM();
+		hmm.readModel(modelFilename);
+
+		
+		
+		PrintStream ps=io.FileUtil.openOutFile(statsFilename);
+		
+		double [][] emit=hmm.getEmitProb();
+		for(int i=0;i<emit.length;i++){
+			ArrayList<IntDoublePair>l=new ArrayList<IntDoublePair>();
+			for(int j=0;j<emit[i].length;j++){
+				l.add(new IntDoublePair(j,emit[i][j]));
+			}
+			Collections.sort(l);
+			ps.println(i);
+			for(int j=0;j<NUM_WORD;j++){
+				if(j>=dict.length){
+					break;
+				}
+				ps.print(dict[l.get(j).idx]+"\t");
+				if((1+j)%10==0){
+					ps.println();
+				}
+			}
+			ps.println("\n");
+		}
+		
+		checkMaxwt(hmm,ps,test.getAllData());
+		
+		int terminalSym=vocab.get(Corpus .END_SYM);
+		//sample 10 sentences
+		for(int i=0;i<10;i++){
+			int []sent=hmm.sample(terminalSym);
+			for(int j=0;j<sent.length;j++){
+				ps.print(dict[sent[j]]+"\t");
+			}
+			ps.println();
+		}
+		
+		ps.close();
+		
+	}
+	
+	public static void checkMaxwt(HMM hmm,PrintStream ps,int [][]data){
+		double [][]emit=hmm.getEmitProb();
+		maxwt=new double[emit.length][emit[0].length];
+		
+		hmm.computeMaxwt(maxwt,data);
+		double sum=0;
+		for(int i=0;i<maxwt.length;i++){
+			for(int j=0;j<maxwt.length;j++){
+				sum+=maxwt[i][j];
+			}
+		}
+		
+		ps.println("max w t P(w_i|t)"+sum);
+		
+	}
+	
+}
diff --git a/gi/posterior-regularisation/prjava/src/test/IntDoublePair.java b/gi/posterior-regularisation/prjava/src/test/IntDoublePair.java
new file mode 100644
index 00000000..3f9f0ad7
--- /dev/null
+++ b/gi/posterior-regularisation/prjava/src/test/IntDoublePair.java
@@ -0,0 +1,23 @@
+package test;
+
+public class IntDoublePair implements Comparable{
+	double val;
+	int idx;
+	public int compareTo(Object o){
+		if(o instanceof IntDoublePair){
+			IntDoublePair pair=(IntDoublePair)o;
+			if(pair.val>val){
+				return 1;
+			}
+			if(pair.val<val){
+				return -1;
+			}
+			return 0;
+		}
+		return -1;
+	}
+	public IntDoublePair(int i,double v){
+		val=v;
+		idx=i;
+	}
+}
diff --git a/gi/posterior-regularisation/prjava/src/test/X2y2WithConstraints.java b/gi/posterior-regularisation/prjava/src/test/X2y2WithConstraints.java
new file mode 100644
index 00000000..9059a59e
--- /dev/null
+++ b/gi/posterior-regularisation/prjava/src/test/X2y2WithConstraints.java
@@ -0,0 +1,131 @@
+package test;
+
+
+
+import optimization.gradientBasedMethods.ProjectedGradientDescent;
+import optimization.gradientBasedMethods.ProjectedObjective;
+import optimization.gradientBasedMethods.stats.OptimizerStats;
+import optimization.linesearch.ArmijoLineSearchMinimizationAlongProjectionArc;
+import optimization.linesearch.InterpolationPickFirstStep;
+import optimization.linesearch.LineSearchMethod;
+import optimization.projections.BoundsProjection;
+import optimization.projections.Projection;
+import optimization.projections.SimplexProjection;
+import optimization.stopCriteria.CompositeStopingCriteria;
+import optimization.stopCriteria.GradientL2Norm;
+import optimization.stopCriteria.ProjectedGradientL2Norm;
+import optimization.stopCriteria.StopingCriteria;
+import optimization.stopCriteria.ValueDifference;
+
+
+/**
+ * @author javg
+ * 
+ * 
+ *ax2+ b(y2 -displacement)
+ */
+public class X2y2WithConstraints extends ProjectedObjective{
+
+
+	double a, b;
+	double dx;
+	double dy;
+	Projection projection;
+	
+	
+	public X2y2WithConstraints(double a, double b, double[] params, double dx, double dy, Projection proj){
+		//projection = new BoundsProjection(0.2,Double.MAX_VALUE);
+		super();
+		projection = proj;	
+		this.a = a;
+		this.b = b;
+		this.dx = dx;
+		this.dy = dy;
+		setInitialParameters(params);
+		System.out.println("Function " +a+"(x-"+dx+")^2 + "+b+"(y-"+dy+")^2");
+		System.out.println("Gradient " +(2*a)+"(x-"+dx+") ; "+(b*2)+"(y-"+dy+")");
+		printParameters();
+		projection.project(parameters);
+		printParameters();
+		gradient = new double[2];
+	}
+	
+	public double getValue() {
+		functionCalls++;
+		return a*(parameters[0]-dx)*(parameters[0]-dx)+b*((parameters[1]-dy)*(parameters[1]-dy));
+	}
+
+	public double[] getGradient() {
+		if(gradient == null){
+			gradient = new double[2];
+		}
+		gradientCalls++;
+		gradient[0]=2*a*(parameters[0]-dx);
+		gradient[1]=2*b*(parameters[1]-dy);
+		return gradient;
+	}
+	
+	
+	public double[] projectPoint(double[] point) {
+		double[] newPoint = point.clone();
+		projection.project(newPoint);
+		return newPoint;
+	}	
+	
+	public void optimizeWithProjectedGradientDescent(LineSearchMethod ls, OptimizerStats stats, X2y2WithConstraints o){
+		ProjectedGradientDescent optimizer = new ProjectedGradientDescent(ls);
+		StopingCriteria stopGrad = new ProjectedGradientL2Norm(0.001);
+		StopingCriteria stopValue = new ValueDifference(0.001);
+		CompositeStopingCriteria compositeStop = new CompositeStopingCriteria();
+		compositeStop.add(stopGrad);
+		compositeStop.add(stopValue);
+		
+		optimizer.setMaxIterations(5);
+		boolean succed = optimizer.optimize(o,stats,compositeStop);
+		System.out.println("Ended optimzation Projected Gradient Descent\n" + stats.prettyPrint(1));
+		System.out.println("Solution: " + " x0 " + o.parameters[0]+ " x1 " + o.parameters[1]);
+		if(succed){
+			System.out.println("Ended optimization in " + optimizer.getCurrentIteration());
+		}else{
+			System.out.println("Failed to optimize");
+		}
+	}
+	
+	
+	
+	public String toString(){
+		
+		return "P1: " + parameters[0] + " P2: " + parameters[1] + " value " + getValue() + " grad (" + getGradient()[0] + ":" + getGradient()[1]+")";
+	}
+	
+	public static void main(String[] args) {
+		double a = 1;
+		double b=1;
+		double x0 = 0;
+		double y0  =1;
+		double dx = 0.5;
+		double dy = 0.2	;
+		double [] parameters = new double[2];
+		parameters[0] = x0;
+		parameters[1] = y0;
+		X2y2WithConstraints o = new X2y2WithConstraints(a,b,parameters,dx,dy, 
+				new SimplexProjection(0.5)
+				//new BoundsProjection(0.0,0.4)
+		);
+		System.out.println("Starting optimization " + " x0 " + o.parameters[0]+ " x1 " + o.parameters[1] + " a " + a + " b "+b );
+		o.setDebugLevel(4);
+		
+		LineSearchMethod ls = new ArmijoLineSearchMinimizationAlongProjectionArc(new InterpolationPickFirstStep(1));
+		
+		OptimizerStats stats = new OptimizerStats();
+		o.optimizeWithProjectedGradientDescent(ls, stats, o);
+		
+//		o = new x2y2WithConstraints(a,b,x0,y0,dx,dy);
+//		stats = new OptimizerStats();
+//		o.optimizeWithSpectralProjectedGradientDescent(stats, o);
+	}
+	
+	
+	
+	
+}