2657 Productions News

I’m not at all religious, but…

Ananda — Sun, 26 Jun 2011 12:43:03 +0000

… here is a goddess that I am happy to worship…

Don't make me squirt my milk bottle at you!

I also had a few alternatives–and I’m still not sure which one is my favorite.

On the trucks around town…

Ananda — Wed, 15 Jun 2011 04:40:22 +0000

Anyone who has spent some time in India is sure to have noticed the slogans painted on the back of trucks, autos, and other vehicles advising “we two, ours one”. This is part of India’s “family planning” efforts–efforts which have had a pretty bumpy history that included a forced sterilization program.

Originally, the slogans were “we two, ours two”, or at least that was the catchy English version–regional languages usually had a slogan more along the lines of “one family, two children”. And, the change to the new slogan led to at least one humorous math discussion with an auto driver who commented that, “Earlier, it was ‘we two, ours two’; now, it is ‘we two, ours one’. What’s next? ‘We two, ours half?’”

Anyway, keen observers might have noticed the following new addition to selected trucks:

Stratified random sampling in R from a data frame

Ananda — Fri, 20 May 2011 18:01:21 +0000

After a little bit more work, there’s a new stratified random sampling function, this one letting you sample from a data frame, returning all the variables for each of your samples as a nice data frame that you can continue working on as usual.

Get the function at http://news.mrdwab.com/stratified. Usage notes in the head of the function.

Here’s the function:

^?View Code RSPLUS

stratified = function(df, id, group, size, seed="NULL", ...) {
  #  USE: * Specify your data frame, ID variable (as column number), and
  #         grouping variable (as column number) as the first three arguments.
  #       * Decide on your sample size. For a sample proportional to the
  #         population, enter "size" as a decimal. For an equal number of
  #         samples from each group, enter "size" as a whole number.
  #       * Decide on if you want to use a seed or not. If not, leave blank
  #         or type "NULL" (with quotes). 
  #
  #  Example 1: To sample 10% of each group from a data frame named "z", where
  #             the ID variable is the first variable, the grouping variable
  #             is the fourth variable, and the desired seed is "1", use:
  # 
  #                 > stratified(z, 1, 4, .1, 1)
  #
  #  Example 2: To run the same sample as above but without a seed, use:
  # 
  #                 > stratified(z, 1, 4, .1)
  #
  #  Example 3: To sample 5 from each group from a data frame named "z", where
  #             the ID variable is the first variable, the grouping variable
  #             is the third variable, and the desired seed is 2, use:
  #
  #                 > stratified(z, 1, 3, 5, 2)
  #
  #  NOTE: Not tested on datasets with LOTS of groups or with HUGE
  #        differences in group sizes. Probably INCREDIBLY inefficient.
 
  k = unstack(data.frame(as.vector(df[id]), as.vector(df[group])))
  l = length(k)
  results = vector("list", l)
 
  if (seed == "NULL" & size < 1) {
      for (i in 1:length(k)) {
        N = k[[i]]
        n = round(length(N)*size)
        results[[i]] = list(sample(N, n, ...))
      }
    } else if (seed == "NULL" & size >= 1) {
      for (i in 1:length(k)) {
        N = k[[i]]
        results[[i]] = list(sample(N, size, ...))
      }
    } else if (size < 1) {
      for (i in 1:length(k)) {
        set.seed(seed)
        N = k[[i]]
        n = round(length(N)*size)
        results[[i]] = list(sample(N, n, ...))
      }
    } else if (size >= 1) {
      for (i in 1:length(k)) {
        set.seed(seed)
        N = k[[i]]
        results[[i]] = list(sample(N, size, ...))
      }
    }
  z = data.frame(c(unlist(results)))
  names(z) = names(df[id])
  w = merge(df, z)
  w[order(w[group]), ]
}

And here are some examples of the function in action:

^?View Code RSPLUS

> source("http://news.mrdwab.com/stratified")
> # Make up some data
> A = 1:100
> B = sample(c("AA", "BB", "CC", "DD", "EE"), 100, replace=T)
> C = rnorm(100)
> D = abs(round(rnorm(100), digits=1))
> E = sample(c("CA", "NY", "TX"), 100, replace=T)
> dat = data.frame(A, B, C, D, E)
> # view the first few rows
> head(dat)
  A  B           C   D  E
1 1 CC -0.07870439 0.6 NY
2 2 CC -0.65048634 0.3 TX
3 3 EE  1.02703616 1.3 NY
4 4 BB -1.08696775 0.4 TX
5 5 CC  0.56741795 0.2 CA
6 6 AA -0.46448941 0.5 TX
> # Sample 10% from each group from variable B, no seed
> stratified(dat, 1, 2, .1)
    A  B           C   D  E
2   6 AA -0.46448941 0.5 TX
7  71 AA  1.98128479 2.1 CA
5  53 BB  1.00539398 0.7 NY
10 97 BB  0.68252675 1.9 NY
1   1 CC -0.07870439 0.6 NY
4  42 CC -2.00256854 0.3 TX
8  76 DD -0.84151459 0.2 NY
9  95 DD -0.47276142 0.3 CA
11 99 DD  1.05173419 2.1 TX
3  10 EE -0.69079473 1.1 TX
6  57 EE -0.38210921 1.5 CA
> # Sample 10% from each group from variable E, seed of 1
> stratified(dat, 1, 5, .1, 1)
    A  B          C   D  E
4  33 AA  1.6105099 0.5 CA
7  48 AA  0.3128274 0.6 CA
9  62 DD  0.4673061 0.0 CA
10 86 EE  0.4047880 1.6 CA
3  28 AA -1.6815553 0.3 NY
5  36 AA  0.3307508 0.3 NY
8  53 BB  1.0053940 0.7 NY
1  21 DD  0.5229282 1.2 TX
2  27 BB  0.8678977 0.7 TX
6  44 DD -0.5790353 0.9 TX
> # You can also be verbose if it helps you remember what you're doing
> stratified(df=dat, id=1, group=5, size=.1, seed=1)
    A  B          C   D  E
4  33 AA  1.6105099 0.5 CA
7  48 AA  0.3128274 0.6 CA
9  62 DD  0.4673061 0.0 CA
10 86 EE  0.4047880 1.6 CA
3  28 AA -1.6815553 0.3 NY
5  36 AA  0.3307508 0.3 NY
8  53 BB  1.0053940 0.7 NY
1  21 DD  0.5229282 1.2 TX
2  27 BB  0.8678977 0.7 TX
6  44 DD -0.5790353 0.9 TX

Stratified Random Sampling in R–A Function in Progress

Ananda — Sun, 15 May 2011 10:16:02 +0000

IMPORTANT: This is here mostly to remind me of how I solved my problem. You should read Stratified random sampling in R from a data frame if you really want to use this function.

I know that sampling is quite complex, and I will admit that I know very little about its complexities. Fortunately, software like R lets you draw simple random samples pretty easily, either either with or without replacement. Unfortunately, I could not find any feature to allow me to do simple stratified random sampling, at least not with the features I was looking for. Fortunately again, with a little bit of experimenting, it can be pretty easy to learn how to write functions in R when a direct solution does not present itself.

This post shares my initial “work-in-progress” on writing an R function for stratified sampling.

The problem…

Here’s the minimum that I was hoping for:

I wanted to be able to draw both a proportional sample (which is more common, for it allows you to make generalizations about the population as a whole) as well as a fixed-size sample (which less common, but it is useful for making comparisons across groups).
I often use a seed when sampling, so I wanted that to be a part of the function.
I wanted the output to be the same as if I were to sample from each group individually.
I was hoping that my output could be stored as a new object that I could then reuse (either a list or a data frame, preferably the latter).

My initial searches directed me to Yihui Xie’s page on stratified sampling using tapply(). However, this option did not satisfy my needs. As far as I could figure, it only allowed me to take a fixed sample size. Also, I wasn’t totally satisfied with the output.

Consider the following. In Yihui Xie’s example, there is a difference between the results one would get if they sampled from each group separately, but using the same seed.

^?View Code RSPLUS

> dat = data.frame(x = 1:15, stratum = gl(3, 5))
> dat
    x stratum
1   1       1
2   2       1
3   3       1
4   4       1
5   5       1
6   6       2
7   7       2
8   8       2
9   9       2
10 10       2
11 11       3
12 12       3
13 13       3
14 14       3
15 15       3
> set.seed(1); tapply(dat$x, dat$stratum, sample, size = 3)
$`1`
[1] 2 5 4
 
$`2`
[1] 10  6  8
 
$`3`
[1] 15 13 12
 
> # Compare with what we get when we sample individually:
> set.seed(1); sample(1:5, 3)
[1] 2 5 4
> set.seed(1); sample(6:10, 3)
[1]  7 10  9
> set.seed(1); sample(11:15, 3)
[1] 12 15 14

I’m sure there’s some sampling theory that explains this, or at least something about how R treats its data, but at the moment, that’s beyond my humble level of expertise.

Stratified sampling, Mr. DWAB style…

The solution I arrived at is to use “unstack()” and a few conditional loops to take the samples.

And, without more rambling, here’s what I came up with.

^?View Code RSPLUS

stratified = function(df, size, seed, dframe=FALSE, ...) {
         # USE: Start with a data frame with your cases in one column and your
         # groups in another column. Decide on if you want to use a seed or not. 
         # If not, seed should be "NO" (with quotes). Decide on if you want your
         # output as a data frame or not; by default, dframe is set to "FALSE".
         # To take a sample proportional to the population size in each group,
         # enter "size" as a decimal. Otherwise, enter size as a whole number.
         #
         # Example 1a: To sample 10% of each group from a data frame named "z"
         # and using a seed of "1", use: > stratified(z, .1, 1)
         # Example 1b: To run the same sample as above but display the result as
         # a data frame, use: > stratified(z, .1, 1, T)
         #
         # Example 2: To sample 10% of each group from a data frame named "z"
         # and using no seed, use: > stratified(z, .1, "NO")
         #
         # Example 3: To sample 5 from each group from a data frame named "z"
         # and using a seed of 30, use: > stratified(z, 5, 30)
         #
         # NOTE: Not recommended for datasets with LOTS of groups or with HUGE
         # differences in group sizes.
  k = unstack(df)
  if (dframe == FALSE) {
    if (seed == "NO" & size < 1) {
      for (i in 1:length(k)) {
        N = k[[i]]
        n = round(length(N)*size)
        pre = structure(list("Group" = names(k[i]), "Population Size" =
                             length(k[[i]]), "Sample Size" = n, Seed = seed,
                             Sample = sample(N, n, ...)), class = "power.htest")
        print(pre)
      }
    } else if (seed == "NO" & size >= 1) {
      for (i in 1:length(k)) {
        N = k[[i]]
        pre = structure(list("Group" = names(k[i]), "Population Size" =
                             length(k[[i]]), "Sample Size" = size, Seed = seed,
                             Sample = sample(N, size, ...)),
                             class = "power.htest")
        print(pre)
      }
    } else if (size < 1) {
      for (i in 1:length(k)) {
        set.seed(seed)
        N = k[[i]]
        n = round(length(N)*size)
        pre = structure(list("Group" = names(k[i]), "Population Size" =
                             length(k[[i]]), "Sample Size" = n, Seed = seed,
                             Sample = sample(N, n, ...)), class = "power.htest")
        print(pre)
      }
    } else if (size >= 1) {
      for (i in 1:length(k)) {
        set.seed(seed)
        N = k[[i]]
        pre = structure(list("Group" = names(k[i]), "Population Size" =
                             length(k[[i]]), "Sample Size" = size, Seed = seed,
                             Sample = sample(N, size, ...)),
                             class = "power.htest")
        print(pre)
      }
    }
  } else {
    if (seed == "NO" & size < 1) {
      for (i in 1:length(k)) {
        N = k[[i]]
        n = round(length(N)*size)
        res = data.frame(names(k[i]), sample(N, n, ...))
        names(res) = c("Group", "Samples")
        print(res)
      }
    } else if (seed == "NO" & size >= 1) {
      for (i in 1:length(k)) {
        N = k[[i]]
        res = data.frame(names(k[i]), sample(N, size, ...))
        names(res) = c("Group", "Samples")
        print(res)
      }
    } else if (size < 1) {
      for (i in 1:length(k)) {
        set.seed(seed)
        N = k[[i]]
        n = round(length(N)*size)
        res = data.frame(names(k[i]), sample(N, n, ...))
        names(res) = c("Group", "Samples")
        print(res)
      }
    } else if (size >= 1) {
      for (i in 1:length(k)) {
        set.seed(seed)
        N = k[[i]]
        res = data.frame(names(k[i]), sample(N, size, ...))
        names(res) = c("Group", "Samples")
        print(res)
      }
    }
  }
}

You can load the function by typing:

^?View Code RSPLUS

1	> source("http://news.mrdwab.com/stratified-beta")

And now, to test it…

Let’s generate some dummy data and see what we can come up with. The function takes the following arguments (in the following order):

df: The source data frame, with the first column being the IDs and the second column being the groups.
size: The sample size you want, either as a percentage (for proportional sampling–expressed as a decimal) or as a whole number.
seed: The seed you want to use. If you don’t want to use a seed, enter “NO”.
dframe: What format you want the output in, either a list or a data frame. Defaults to a list (dframe=FALSE), which is better at the moment since the data frame option is not working the way I expect it to yet.

^?View Code RSPLUS

> # Generate some data
> a = 1:100
> set.seed(123)
> b = sample(c("a", "b", "c", "d"), 100, replace = T)
> z = data.frame(a, b)
> # Check how big each group is
> table(z$b)
 
 a  b  c  d
26 27 20 27
> # Make sure the function is loaded before you continue!
> # source("http://news.mrdwab.com/stratified-beta")
> # Take a 15% sample and use a seed of 1
> stratified(z, .15, 1)
 
 
 
          Group = a
Population Size = 26
    Sample Size = 4
           Seed = 1
         Sample = 38, 45, 54, 81
 
 
 
 
          Group = b
Population Size = 27
    Sample Size = 4
           Seed = 1
         Sample = 39, 43, 60, 79
 
 
 
 
          Group = c
Population Size = 20
    Sample Size = 3
           Seed = 1
         Sample = 23, 26, 33
 
 
 
 
          Group = d
Population Size = 27
    Sample Size = 4
           Seed = 1
         Sample = 21, 31, 53, 71
 
> # Take a sample of 5 from each group and use a seed of 1
> stratified(z, 5, 1)
 
 
 
          Group = a
Population Size = 26
    Sample Size = 5
           Seed = 1
         Sample = 38, 45, 54, 81, 30
 
 
 
 
          Group = b
Population Size = 27
    Sample Size = 5
           Seed = 1
         Sample = 39, 43, 60, 79, 19
 
 
 
 
          Group = c
Population Size = 20
    Sample Size = 5
           Seed = 1
         Sample = 23, 26, 33, 78, 14
 
 
 
 
          Group = d
Population Size = 27
    Sample Size = 5
           Seed = 1
         Sample = 21, 31, 53, 71, 11
 
> # Take a sample of 15 from each group, with replacement, and a seed of 1
> stratified(z, 15, 1, replace=T)
 
 
 
          Group = a
Population Size = 26
    Sample Size = 15
           Seed = 1
         Sample = 38, 45, 56, 91, 35, 91, 96, 74, 62, 15, 35, 30, 74, 45, 81
 
 
 
 
          Group = b
Population Size = 27
    Sample Size = 15
           Seed = 1
         Sample = 39, 44, 63, 93, 29, 93, 95, 66, 64, 3, 29, 19, 70, 44, 77
 
 
 
 
          Group = c
Population Size = 20
    Sample Size = 15
           Seed = 1
         Sample = 23, 26, 55, 94, 22, 92, 94, 72, 61, 9, 22, 14, 72, 26, 78
 
 
 
 
          Group = d
Population Size = 27
    Sample Size = 15
           Seed = 1
         Sample = 21, 32, 58, 88, 16, 88, 89, 65, 59, 4, 16, 11, 67, 32, 69
 
> # Take a sample of 10% from each group, using a seed of 1,
> # and display the output as a data frame
> stratified(z, .1, 1, dframe=T)
  Group Samples
1     a      38
2     a      45
3     a      54
  Group Samples
1     b      39
2     b      43
3     b      60
  Group Samples
1     c      23
2     c      26
  Group Samples
1     d      21
2     d      31
3     d      53

Replicating the results from tapply()

I mentioned earlier that the results are different from what you would get if you were to use the tapply() function. However, it is easy to get the same results using this stratified function–simply move your “seed” outside of the function (enter seed as "NO" [with quotes] and instead, use set.seed() as you normally would).

^?View Code RSPLUS

> # See what tapply() gives us
> set.seed(1); tapply(z$a, z$b, sample, size = 4)
$a
[1] 38 45 54 81
 
$b
[1] 29 86 95 63
 
$c
[1] 61  9 14 92
 
$d
[1] 67 31 68 34
 
> # The normal usage for the stratified function
> stratified(z, 4, 1, dframe = T)
  Group Samples
1     a      38
2     a      45
3     a      54
4     a      81
  Group Samples
1     b      39
2     b      43
3     b      60
4     b      79
  Group Samples
1     c      23
2     c      26
3     c      33
4     c      78
  Group Samples
1     d      21
2     d      31
3     d      53
4     d      71
> # Getting the same results as tapply()
> # Set the seed before using the function,
> # and set the seed for the function as "NO"
> set.seed(1); stratified(z, 4, "NO", dframe = T)
  Group Samples
1     a      38
2     a      45
3     a      54
4     a      81
  Group Samples
1     b      29
2     b      86
3     b      95
4     b      63
  Group Samples
1     c      61
2     c       9
3     c      14
4     c      92
  Group Samples
1     d      67
2     d      31
3     d      68
4     d      34

The unfortunate…

There are some advantages to each of the output formats. I’ve set up the list to be quite verbose, which is useful with the proportionate sampling since it shows us how many samples have been taken from each group. The data frame output format, on the other hand, is quite compact.

What I still need to figure out, though, is why R won’t store my output. I suspect that it has something to do with how my loops are set up. I assume that somewhere, I need to add something like an rbind command.

When the time is right, I will be sure to post what I’ve found.

Reshaping data in R revisited

Ananda — Mon, 18 Apr 2011 04:15:17 +0000

A year ago, I wrote a post about reshaping data from a wide format to a long format. I thought that considering how much time had passed, it would be good to revisit R’s in-built reshape functions.

For these examples, I’ve copied the Stata examples from the UCLA Academic Technology Services’s “Reshape data wide to long” page. Since the data is provided in Stata dta files, you need to first load the “foreign” package to be able to read the data in R.

This first example is very basic. There are four variables, the first one being the unique id, and the remaining three being the measures over three years.

The basic reshape command in R needs you to specify the data that is being reshaped, the ultimate “direction” (wide or long), and which variables are the ones to be reshaped. The default character R expects for “sep” is a period–in other words, it expects that your variables are named in the form of “faminc.96″ and so on. However, your variables may be named in other ways, for example “faminc-96″, “faminc_96″, or (as in this example) “faminc96″. If your variable naming pattern is anything other than what R expects as its default, you also need to specify the separating character. In the case of this dataset, there is no separating character, so you simply use sep="".

^?View Code RSPLUS

> library(foreign) # Lets us use Stata files directly
> # Example 1: Very basic reshape
> # Use "read.dta" instead of "read.csv" or "read.table"
> faminc = read.dta("http://www.ats.ucla.edu/stat/stata/modules/faminc.dta")
> faminc
  famid faminc96 faminc97 faminc98
1     3    75000    76000    77000
2     1    40000    40500    41000
3     2    45000    45400    45800
> l.faminc = reshape(faminc, direction="long", varying=2:4, sep="", idvar="famid")
> l.faminc[order(l.faminc$famid),]
     famid time faminc
1.96     1   96  40000
1.97     1   97  40500
1.98     1   98  41000
2.96     2   96  45000
2.97     2   97  45400
2.98     2   98  45800
3.96     3   96  75000
3.97     3   97  76000
3.98     3   98  77000

In the second example at UCLA ATS, the unique identifier is the combination of the first two variables. Since R assumes that whatever you have not specified as varying is going to be your identifying variable, it is not always required that you specify anything for “idvar”.

^?View Code RSPLUS

> # Example 2: Two identifying variables
> kidshtwt = read.dta("http://www.ats.ucla.edu/stat/stata/modules/kidshtwt.dta")
> kidshtwt
  famid birth ht1 ht2 wt1 wt2
1     1     1 2.8 3.4  19  28
2     1     2 2.9 3.8  21  28
3     1     3 2.2 2.9  20  23
4     2     1 2.0 3.2  25  30
5     2     2 1.8 2.8  20  33
6     2     3 1.9 2.4  22  33
7     3     1 2.2 3.3  22  28
8     3     2 2.3 3.4  20  30
9     3     3 2.1 2.9  22  31
> # Note the use of "timevar" to name the "times" column more appropriately.
> # Also, we want to exclude the weight data from our reshape.
> l.kidsht = reshape(kidshtwt[-c(5, 6)], direction="long", varying=3:4, sep="",
+                    idvar=1:2, timevar="age")
> # There are other ways to do that previous step. Here is one:
> #      > l.kidsht = reshape(kidshtwt, direction="long", idvar=1:2,
> #      +                    varying=3:4, drop=5:6, sep="", timevar="age")
> # Let's sort the data first by the family id then by the birth order
> l.kidsht[order(l.kidsht$famid, l.kidsht$birth),]
      famid birth age  ht
1.1.1     1     1   1 2.8
1.1.2     1     1   2 3.4
1.2.1     1     2   1 2.9
1.2.2     1     2   2 3.8
1.3.1     1     3   1 2.2
1.3.2     1     3   2 2.9
2.1.1     2     1   1 2.0
2.1.2     2     1   2 3.2
2.2.1     2     2   1 1.8
2.2.2     2     2   2 2.8
2.3.1     2     3   1 1.9
2.3.2     2     3   2 2.4
3.1.1     3     1   1 2.2
3.1.2     3     1   2 3.3
3.2.1     3     2   1 2.3
3.2.2     3     2   2 3.4
3.3.1     3     3   1 2.1
3.3.2     3     3   2 2.9

The third example at UCLA’s page is pretty straightforward. It uses the same data we just loaded in example 2, but we are reshaping all four measured variables.

^?View Code RSPLUS

> l.kidshtwt = reshape(kidshtwt, direction="long", idvar=1:2,varying=3:6,
+                      sep="", timevar="age")
> l.kidshtwt[order(l.kidshtwt$famid, l.kidshtwt$birth),]
      famid birth age  ht wt
1.1.1     1     1   1 2.8 19
1.1.2     1     1   2 3.4 28
1.2.1     1     2   1 2.9 21
1.2.2     1     2   2 3.8 28
1.3.1     1     3   1 2.2 20
1.3.2     1     3   2 2.9 23
2.1.1     2     1   1 2.0 25
2.1.2     2     1   2 3.2 30
2.2.1     2     2   1 1.8 20
2.2.2     2     2   2 2.8 33
2.3.1     2     3   1 1.9 22
2.3.2     2     3   2 2.4 33
3.1.1     3     1   1 2.2 22
3.1.2     3     1   2 3.3 28
3.2.1     3     2   1 2.3 20
3.2.2     3     2   2 3.4 30
3.3.1     3     3   1 2.1 22
3.3.2     3     3   2 2.9 31

The fourth example was the most tricky one for me at first. In that example, the variables are not distinuished by “time” (numerically) but rather, by a character. As you can see, the variable names are “famid”, “named”, “incd”, “namem”, and “incm” — in other words, income and name for dad (variable ending in “d”) and mom (variable ending in “m”) for each family.

^?View Code RSPLUS

> # Example 3: Non-numeric identifiers for the variables
> dadmomw = read.dta("http://www.ats.ucla.edu/stat/stata/modules/dadmomw.dta")
> dadmomw
  famid named  incd namem  incm
1     1  Bill 30000  Bess 15000
2     2   Art 22000   Amy 18000
3     3  Paul 25000   Pat 50000
> r.dadmomw = reshape(dadmomw, direction="long", idvar=1, varying=2:5,
+                     sep="", v.names=c("name", "inc"),
+                     timevar="dadmom", times=c("dad", "mom"))
> r.dadmomw[order(r.dadmomw$famid),]
      famid dadmom  name  inc
1.dad     1    dad 30000 Bill
1.mom     1    mom 15000 Bess
2.dad     2    dad 22000  Art
2.mom     2    mom 18000  Amy
3.dad     3    dad 25000 Paul
3.mom     3    mom 50000  Pat

Stata’s commands are certainly more direct (see below for what you would do for the last example in Stata). R’s commands sometimes tend to be a bit verbose, but in some ways, that might also help you remember what you’re doing. (I still don’t know what “i” and “j” in the Stata reshape commands stand for.) If you can afford the ~ $2,500 price tag, Stata is also faster.

^?View Code STATA

use http://www.ats.ucla.edu/stat/stata/modules/dadmomw, clear
list
use http://www.ats.ucla.edu/stat/stata/modules/dadmomw, clear
list
reshape long name  inc, i(famid) j(dadmom) string
list