Q1. Git/GitHub

No handwritten homework reports are accepted for this course. We work with Git and GitHub. Efficient and abundant use of Git, e.g., frequent and well-documented commits, is an important criterion for grading your homework.

  1. Apply for the Student Developer Pack at GitHub using your UCLA email.

  2. Create a private repository biostat-m280-2018-winter and add Hua-Zhou and juhkim111 as your collaborators with write permission.

  3. Top directories of the repository should be hw1, hw2, … Maintain two branches master and develop. The develop branch will be your main playground, the place where you develop solution (code) to homework problems and write up report. The master branch will be your presentation area. Submit your homework files (R markdown file Rmd, html file converted from R markdown, all code and data sets to reproduce results) in master branch.

  4. After each homework due date, teaching assistant and instructor will check out your master branch for grading. Tag each of your homework submissions with tag names hw1, hw2, … Tagging time will be used as your submission time. That means if you tag your hw1 submission after deadline, penalty points will be deducted for late submission.

Q2. Linux Shell Commands

The 35.227.165.60:/home/m280-data/hw1 folder contains a typical genetic data set in plink format. If interested, you can read plink documentation at http://zzz.bwh.harvard.edu/plink/. But it’s definitely not necessary for this homework.

Please, do not put these data files into Git; they are huge. You even don’t need to copy them into your directory. Just read from the data folder /home/m280-data/hw1 directly.

Use Linux shell commands to answer following questions.

  1. How many persons are in the data set (statisticians call this n)? How many SNPs are in the data set (statisticians call this p)?

  2. Which chromosomes does this data set contain? How many SNPs are in each chromosome?

  3. MAP4 (microtubule-associated protein 4) is a gene on chromosome 3 spanning positions 47,892,180 bp – 48,130,769 bp. How many SNPs are located within MAP4 gene?

  4. Statistical geneticists often have to reformat a data set to feed into various analysis programs. For example, to use the Mendel software http://www.genetics.ucla.edu/software/mendel, we have to reformat the data set to be read by Mendel.

    • Mendel’s SNP definition file is similar to the plink bim file but has format
      SNP ID, Chromosome, Base Pair Position
      with each field separated by a comma. Write a Linux shell command to convert merge-geno.bim to Mendel SNP definition file. The first few lines of the Mendel SNP definition file should look like
    ##     2.40 = FILE FORMAT VERSION NUMBER.
## 8348674  = NUMBER OF SNPS LISTED HERE.
## 1-54490,1,54490
## 1-55550,1,55550
## 1-57033,1,57033
## 1-57064,1,57064
## 1-57818,1,57818
## 1-58432,1,58432
## 1-58448,1,58448
## 1-58814,1,58814
    • Mendel’s pedigree file is similar to the plink fam file but has format
      Family ID, Person ID, Father ID, Mother ID, Sex coded as M or F, Twin Status
      with each field separated by a comma. Write a Linux shell command to convert merge-geno.fam to Mendel pedigree file. Since twin status is not available in plink format, we put nothing for that field. Also Mendel limits Person ID to have length less than or equal to 8 characters, so we have to strip the string T2DG from the IDs. First few lines of the Mendel pedigree should look like
    ## 2,0200001,,,M,
## 2,0200002,,,F,
## 2,0200003,,,F,
## 2,0200004,,,F,
## 2,0200005,,,M,
## 2,0200006,,,M,
## 2,0200007,,,F,
## 2,0200008,,,F,
## 2,0200009,,,F,
## 2,0200012,,,M,
## 2,0200013,,,M,
## 2,0200018,,,M,
## 2,0200023,,,F,
## 2,0200024,,,M,
## 2,0200027,,,F,
## 2,0200031,0200001,0200015,M,
## 2,0200032,0200001,0200015,F,
## 2,0200033,0200001,0200015,F,
## 2,0200034,0200001,0200015,F,
## 2,0200035,0200001,0200015,F,

Q3. R Batch Run

In class we discussed using R to organize simulation studies.

  1. Expand the runSim.R script to include arguments seed (random seed), n (sample size), dist (distribution) and rep (number of simulation replicates). When dist="gaussian", generate data from standard normal; when dist="t1", generate data from t-distribution with degree of freedom 1 (same as Cauchy distribution); when dist="t5", generate data from t-distribution with degree of freedom 5. Calling runSim.R will (1) set random seed according to argument seed, (2) generate data according to argument dist, (3) compute the primed-indexed average estimator in class and the classical sample average estimator for each simulation replicate, (4) report the average mean squared error (MSE) \[ \frac{\sum_{r=1}^{\text{rep}} (\widehat \mu_r - \mu_{\text{true}})^2}{\text{rep}} \] for both methods.

  2. Modify the autoSim.R script to run simulations with combinations of sample sizes nVals = seq(100, 500, by=100) and distributions distTypes = c("gaussian", "t1", "t5") and write output to appropriately named files. Use rep = 50, and seed = 280.

  3. Write an R script to collect simulation results from output files and print average MSEs in a table of format

\(n\) Method \(t_1\) \(t_5\) Gaussian
100 PrimeAvg
SampAvg
200 PrimeAvg
SampAvg
300 PrimeAvg
SampAvg
400 PrimeAvg
SampAvg
500 PrimeAvg
SampAvg