Integration of sampling weights in a factor model

4 replies [Last post]
metavid's picture
Offline
Joined: 06/23/2014

Good day,

I saw a thread discussing about integration of sample weights with likelihood but I still was not able to work out the code given in that thread. The code I used and worked out is stated below (from an earlier thread):

fullModel <- mxModel("ThisIsHowYouDoWeights",
mxModel("firstModel",
....
mxFIMLObjective("myCov", "myMeans", dims, vector=TRUE)),
mxData(myData, "raw"),
mxAlgebra(-2 * sum(data.weight %x% log(firstModel.objective), name="obj"),
mxAlgebraObjective("obj")
)

with data.weight being a vector (in my case) and the firstModel.objective (also a vector) is calculated differently (not with vector=TRUE) because it was not working for me.

My analysis is about estimation of parameters from a one-factor model (with sampling weights). The model without the sampling weights (normal model) works alright but there were some issues with extracting the likelihood vector when vector=TRUE on mxFIMLObjective. I then decided to use another code from the website of OpenMx (using the formula of row-likelihood) and I was able to extract the likelihoods with it. Now, with available sample weights and the vector of likelihoods, I was wondering if there is another way to tackle this issue? Also, in my case, because data.weight is a vector as well the firstModel.objective, I multiplied them both (*) instead of using Kronecker multiplication (I thought it was correct?).

Thank you in advance for your suggestions and tips.

metavid

tbates's picture
Offline
Joined: 07/31/2009
need to clarify why vector=TRUE not used in the unweighted model

To generate a weighted likelihood, you need row likelihoods, and to weight and sum each of these to form a single number (likelihood) that can be the ultimate objective of.

You say that firstModel.objective "is calculated differently (not with vector=TRUE) because it was not working for me."

You have to get a row-vector out of that first model. It is those individual likelihoods for each row that will be weighted. So you need to get that working.
Post here for help (post the full model script, and the data if you can).

fullModel <- mxModel("ThisIsHowYouDoWeights",
   mxModel("firstModel", ...
      mxFIMLObjective("myCov", "myMeans", dims, vector=TRUE)
   ),
   mxData(myData, "raw"),
   mxAlgebra(-2 * sum(data.weight %x% log(firstModel.objective), name = "obj"),
   mxAlgebraObjective("obj")
)

You also need to ensure that "myData" contains the column "weight"

> The unweighted model works, but there were issues with extracting the likelihood vector when vector=TRUE

that's the problem to be clear about: what do you mean "issues"? firstModel.objective didn't contain the row likelihoods?

> I then decided to use another code from the website
> (using the formula of row-likelihood)
> and I was able to extract the likelihoods with it.

Which was that?

metavid's picture
Offline
Joined: 06/23/2014
Clarification and questions

Good day,
Firstly, Thank you for providing a response.

I attached a portion of the data set (samp.txt), the data set is much longer but even with this smaller version, I was not able to extract the likelihood vector.

At the moment, the biggest issue is getting the likelihood vector when vector=TRUE. In my previous analysis, I did not incorporate the weights in myData so I am going to do that in the renewed analysis. The labels of the parameters are not fully specified in the code but I don't think that is the reason of the issue.

The code is:

library(foreign)
samp <- read.table("samp.txt",header=TRUE)
nsamp <- names(samp)
example1 <- mxModel("examp1",
# factor loadings
mxMatrix("Full", nrow=6, ncol=1,
values=c(1,rep(0.8,5)),
free=c(FALSE,rep(TRUE,5)), lbound = c(NA,rep(0.01,5)),
name="A"),
# specify the factor intercorrelation matrix
mxMatrix("Symm", nrow=1, ncol=1, values=2, free=TRUE, lbound = 0.01, name="L"),
# specify the matrix of unique item variances
mxMatrix("Diag", nrow=6, ncol=6, values=rep(0.8,6), free=TRUE, lbound = rep(0.01,6), name="U"),
# specify the algebra that results in the model expectations
mxAlgebra(A %*% L %*% t(A) + U,
name="Covar1"
),
# specify a model for the means
mxMatrix("Full", nrow=1, ncol=6,
values=c(rep(3.5,6)), free=TRUE, lbound = rep(0.01,6),
name="Mean1"
),
# choose the full information maximum likelihood objective function
mxFIMLObjective(covariance="Covar1", means="Mean1",dimnames=nsamp,vector=TRUE),
# attach the data to the model
mxData(samp, type="raw")
)

model1 <- mxRun(example1)
stat <- summary(model)

It then gives a warning message:
In model 'examp1' NPSOL returned a non-zero status code 1. The final iterate satisfies the optimality conditions to the accuracy requested, but the sequence of iterates has not yet converged. NPSOL was terminated because no further improvement could be made in the merit function (Mx status GREEN).

I saw in a number of threads that there are a number of factors as to why the warning measage emerges, e.g. modifying starting values and providing lower bounds. I tried to use differing starting values and lower bounds but the warning message is still present.

As to getting the likelihood vector from another source, I used the code from this OpenMX link about FIML specification (I did not put the code because it could be too long): http://openmx.psyc.virginia.edu/docs/openmx/latest/fiml_rowobjective.html. It is underneath the page, starting with model specification and ending with model fitting. I extracted the likelihood vector from one of its output.

I hope I provided sufficient information to tackle this problem. Thank you.

metavid

AttachmentSize
samp.txt 2.68 KB
mhunter's picture
Offline
Joined: 07/31/2009
Example

Here is a full working example

samp <- read.table("samp.txt",header=TRUE)
nsamp <- names(samp)
samp$weight <- runif(nrow(samp)) # Your sample weights!!!
 
example1 <- mxModel("examp1",
	# factor loadings
	mxMatrix("Full", nrow=6, ncol=1,
		values=c(1,rep(0.8,5)),
		free=c(FALSE,rep(TRUE,5)), lbound = c(NA,rep(0.01,5)),
		name="A"),
	# specify the factor intercorrelation matrix
	mxMatrix("Symm", nrow=1, ncol=1, values=2, free=TRUE, lbound = 0.01, name="L"),
	# specify the matrix of unique item variances
	mxMatrix("Diag", nrow=6, ncol=6, values=rep(0.8,6), free=TRUE, lbound = rep(0.01,6), name="U"),
	# specify the algebra that results in the model expectations
	mxAlgebra(A %*% L %*% t(A) + U,
		name="Covar1"
	),
	# specify a model for the means
	mxMatrix("Full", nrow=1, ncol=6,
		values=c(rep(3.5,6)), free=TRUE, lbound = rep(0.01,6),
		name="Mean1"
	),
	# choose the full information maximum likelihood objective function
	mxFIMLObjective(covariance="Covar1", means="Mean1",dimnames=nsamp,vector=TRUE)
)
 
example1a <- mxModel(model="Container", example1,
	# attach the data to the model
	mxData(samp, type="raw"),
	# Algebra for weights
	mxAlgebra(-2 * sum(data.weight %x% log(examp1.objective)), name = "obj"),
	mxAlgebraObjective("obj")
)
 
model1a <- mxRun(example1a)
stat <- summary(model1a)
stat

metavid's picture
Offline
Joined: 06/23/2014
Solved

It works. Again, thank you for providing a solution to my second question!