#88 moving MVNORM models into subfolder. Added MixedModel

.gitignore

@@ -120,3 +120,4 @@ demo/temp.init.R
 ddToolbox/models/parametric_models/exponential_power_models/separateExpPower
 
 ddToolbox/models/parametric_models/exponential_power_models/mixedExpPower
+demo/jags.dump1.R

ddToolbox/models/parametric_models/hyperbolic_magnitude_effect_models/MVNORM_models/ModelMixedME_MVNORM.m

@@ -0,0 +1,29 @@
+classdef ModelMixedME_MVNORM < Hyperbolic1MagEffectMV
+	%ModelHierarchicalME A model to estimate the magnitide effect
+	%   Detailed explanation goes here
+
+	methods (Access = public, Hidden = true)
+
+		function obj = ModelMixedME_MVNORM(data, varargin)
+			obj = obj@Hyperbolic1MagEffectMV(data, varargin{:});
+			obj.modelFilename = 'mixedMEmvnorm';
+			obj = obj.addUnobservedParticipant('GROUP');
+
+			% MUST CALL THIS METHOD AT THE END OF ALL MODEL-SUBCLASS CONSTRUCTORS
+			obj = obj.conductInference();
+		end
+
+    end
+
+    methods (Access = protected)
+
+		function initialParams = initialiseChainValues(obj, nchains) 
+			% Generate initial values of the root nodes
+			for chain = 1:nchains
+				initialParams(chain).groupW			= rand;
+				initialParams(chain).groupALPHAmu	= rand*10;
+				initialParams(chain).groupALPHAsigma= rand*10;
+			end
+		end
+	end
+end

ddToolbox/models/parametric_models/hyperbolic_magnitude_effect_models/MVNORM_models/mixedMEmvnorm.jags

@@ -0,0 +1,82 @@
+
+model{
+
+# DISCOUNT FUNCTION PARAMETERS =================================================
+# RANDOM (BY PARTICIPANT) FACTORS; HYPER-PRIORS = YES
+
+# hyperpriors for (m,c)
+
+
+for (p in 1:(nRealExperimentFiles+1)){ # +1 for unobserved participant
+    mc_mu[p,1] ~ dnorm(-0.243, 1/( (0.027)^2)) # prior over m_mu
+    mc_mu[p,2] ~ dnorm(-3    , 1/( 2^2) )      # prior over c_mu
+    precision[p,1] ~ dgamma(0.001,0.001)
+    precision[p,2] ~ dgamma(0.001,0.001)
+    r[p] ~ dunif(-1,1)
+
+    # reparameterise
+    sigma[p,1] <- 1/sqrt(precision[p,1])
+    sigma[p,2] <- 1/sqrt(precision[p,2])
+    T[p,1,1] <- 1/precision[p,1]
+    T[p,1,2] <- r[p]*sigma[p,1]*sigma[p,2]
+    T[p,2,1] <- r[p]*sigma[p,1]*sigma[p,2]
+    T[p,2,2] <- 1/precision[p,2]
+    TI[p,1:2,1:2] <- inverse(T[p,1:2,1:2])
+
+    mc[p,1:2] ~ dmnorm(mc_mu[p,1:2],TI[p,1:2,1:2])
+
+    m[p] <- mc[p,1]
+    c[p] <- mc[p,2]
+}
+
+# MODEL-SPECIFIC: CALCULATION OF PRESENT SUBJECTIVE VALUES
+for (t in 1:length(ID)) {
+	# Calculate log discount rate for each reward
+	lkA[t] <- m[ID[t]]*log(abs(A[t]))+c[ID[t]]
+	lkB[t] <- m[ID[t]]*log(abs(B[t]))+c[ID[t]]
+
+	# calculate present subjective value for each reward
+	VA[t] <- A[t] / (1+(exp(lkA[t])*DA[t]))
+	VB[t] <- B[t] / (1+(exp(lkB[t])*DB[t]))
+}
+
+# RESPONSE ERROR PARAMETERS ====================================================
+# comparison acuity (alpha)
+groupALPHAmu        ~ dunif(0,1000)
+groupALPHAsigma     ~ dunif(0,1000)
+
+# error rates (epsilon)
+groupW          ~ dbeta(1.1, 10.9)  # mode for lapse rate
+groupKminus2    ~ dgamma(0.01,0.01) # concentration parameter
+groupK          <- groupKminus2+2
+
+epsilon_alpha <- groupW*(groupK-2)+1
+epsilon_beta <- (1-groupW)*(groupK-2)+1
+
+for (p in 1:(nRealExperimentFiles+1)){ # +1 for unobserved participant
+	epsilon[p] ~ dbeta(epsilon_alpha , epsilon_beta ) T(,0.5)
+
+    # using reparameterisation to avoid funnel of hell
+    alpha_offset[p] ~ dnorm(0,1) T(0,)
+    alpha[p] <- groupALPHAmu + alpha_offset[p] * groupALPHAsigma
+}
+
+# MODEL IN-SPECIFIC CODE BELOW... SHOULD NOT CHANGE ACROSS MODELS ==============
+
+# Psychometric function
+for (t in 1:length(ID)) {
+	P[t] <- epsilon[ID[t]] + (1-2*epsilon[ID[t]]) * phi( (VB[t]-VA[t]) / alpha[ID[t]] )
+}
+
+# response likelihood
+for (t in 1:length(ID)) {
+	R[t]  ~ dbern(P[t]) # likelihood of actual response
+    log_lik[t] <- logdensity.bern(R[t], P[t])
+}
+
+# POSTERIOR PREDICTION
+for (t in 1:length(ID)) {
+	Rpostpred[t] ~ dbern(P[t])
+}
+
+}