How to test if two populations have the same demographics? I have 2 population samples for which I want to assess whether they have comparable demographics. What is the right way to test this, that includes multivariate combinations ?
I.e I would like to check that both datasets have:


*

*the same proportion of men

*the same proportion of age [0-18]

*the same proportion of men aged [0-18]
etc.


Here is an example of the data. Ideally, I would like to do the test in R.

 A: Considering that you have 4 categorical variables, you could run a mosaic plot. There is a mosaic function in the vcd library. So, something like the following would give a simple plot (untested code):
install.packages("vcd")
library(vcd)
mosaic(dataframename, main = "Title")

It is relatively simple to add e.g. Pearson residuals, see the help page for mosaic. 
Alternatively, you could do logistic regression with "group" as the dependent variable and the demographic variables as independent variables. If you like, you could then include interactions. 
A: Here is a complete solution in R based on Peter Flom's second suggestion:
First a synthetic dataset:
data <- structure(list(age = structure(c(1L, 2L, 3L, 4L, 1L, 2L, 3L, 
4L, 1L, 2L, 3L, 4L, 1L, 2L, 3L, 4L), .Label = c("0-18", "19-35", 
"36-55", "56+"), class = "factor"), gender = structure(c(2L, 
2L, 2L, 2L, 1L, 1L, 1L, 1L, 2L, 2L, 2L, 2L, 1L, 1L, 1L, 1L), .Label = c("female", 
"male"), class = "factor"), region = structure(c(1L, 1L, 1L, 
1L, 1L, 1L, 1L, 1L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L), .Label = c("A", 
"B"), class = "factor"), n_sample1 = c(546, 520, 409, 461, 458, 
460, 507, 447, 482, 588, 514, 416, 542, 595, 432, 476), n_sample2 = c(186, 
174, 146, 155, 155, 161, 169, 153, 163, 203, 171, 146, 182, 205, 
149, 162)), .Names = c("age", "gender", "region", "n_sample1", 
"n_sample2"), class = "data.frame", row.names = c(NA, -16L))

Now the logistic regressions:
# first order terms
glm.res <- glm(cbind(n_sample2, n_sample1)~age+gender+region,
               data = data, family = "binomial")
summary(glm.res)

# first + second order
glm.res <- glm(cbind(n_sample2, n_sample1)~(age+gender+region)^2,
               data = data, family = "binomial")
summary(glm.res)

# first, second and third order
glm.res <- glm(cbind(n_sample2, n_sample1)~(age+gender+region)^3,
               data = data, family = "binomial")
summary(glm.res)

