<commit_before>package org.dizhang.seqspark.assoc

import breeze.linalg._
import org.dizhang.seqspark.stat.ScoreTest.NullModel
import org.dizhang.seqspark.stat.{MultivariateNormal, Resampling, ScoreTest}
import org.dizhang.seqspark.util.General.RichDouble

import scala.language.existentials
import scala.util.{Success, Try}

/**
  * Variable threshold association method
  */
@SerialVersionUID(7727880001L)
trait VT extends AssocMethod {
  def nullModel: NullModel
  def x: Encode.VT
  def result: AssocMethod.Result
}

object VT {

  def apply(nullModel: NullModel,
            x: Encode.Coding): AnalyticTest = {
    AnalyticTest(nullModel, x.asInstanceOf[Encode.VT])
  }

  def apply(ref: Double, min: Int, max: Int,
            nullModel: NullModel,
            x: Encode.Coding): ResamplingTest = {
    ResamplingTest(ref, min, max, nullModel, x.asInstanceOf[Encode.VT])
  }

  def getStatistic(nm: NullModel, x: Encode.Coding): Double = {
    //println(s"scores: ${st.score.toArray.mkString(",")}")
    //println(s"variances: ${diag(st.variance).toArray.mkString(",")}")
    val m = x.asInstanceOf[Encode.VT].coding
    val ts = m.map{sv =>
      val st = ScoreTest(nm, sv)
      st.score(0)/st.variance(0, 0).sqrt
    }
    //val ts = st.score :/ diag(st.variance).map(x => x.sqrt)
    max(ts)
  }

  @SerialVersionUID(7727880101L)
  final case class AnalyticTest(nullModel: NullModel,
                                x: Encode.VT)
    extends VT with AssocMethod.AnalyticTest
  {
    val scoreTest = {
      //val cnt = x.coding.activeSize
      //println(s"geno: ${geno.coding.rows} x ${geno.coding.cols}, " +
      //  s"active: $cnt sample: ${geno.coding.toDense(0,::).t.toArray.mkString(",")}")
      val m = DenseVector.horzcat(x.coding:_*)
      ScoreTest(nullModel, m)
    }
    val statistic = getStatistic(nullModel, x)
    private val ss = diag(scoreTest.variance).map(_.sqrt)
    private val ts = scoreTest.score :/ ss
    val pValue = {
      val maxT = max(ts)
      val cutoff = maxT * ss
      val pTry = Try(MultivariateNormal.Centered(scoreTest.variance).cdf(cutoff).pvalue)
      pTry match {
        case Success(p) => Some(1.0 - p)
        case _ => None
      }
    }
    def result: AssocMethod.VTAnalytic = {
      val info = s"Ts=${ts.toArray.mkString(",")}"
      AssocMethod.VTAnalytic(x.vars, x.size, statistic, pValue, info)
    }
  }

  @SerialVersionUID(7727880201L)
  final case class ResamplingTest(refStatistic: Double,
                                  min: Int,
                                  max: Int,
                                  nullModel: NullModel,
                                  x: Encode.VT)
    extends VT with AssocMethod.ResamplingTest
  {
    def pCount = {
      Resampling.Simple(refStatistic, min, max, nullModel, x, getStatistic).pCount
    }
    def result: AssocMethod.VTResampling =
      AssocMethod.VTResampling(x.vars, x.size, refStatistic, pCount)
  }

}
<commit_msg>fix bug in analytic VT<commit_after>package org.dizhang.seqspark.assoc

import breeze.linalg._
import org.dizhang.seqspark.stat.ScoreTest.NullModel
import org.dizhang.seqspark.stat.{MultivariateNormal, Resampling, ScoreTest}
import org.dizhang.seqspark.util.General.RichDouble

import scala.language.existentials
import scala.util.{Success, Try}

/**
  * Variable threshold association method
  */
@SerialVersionUID(7727880001L)
trait VT extends AssocMethod {
  def nullModel: NullModel
  def x: Encode.VT
  def result: AssocMethod.Result
}

object VT {

  def apply(nullModel: NullModel,
            x: Encode.Coding): AnalyticTest = {
    AnalyticTest(nullModel, x.asInstanceOf[Encode.VT])
  }

  def apply(ref: Double, min: Int, max: Int,
            nullModel: NullModel,
            x: Encode.Coding): ResamplingTest = {
    ResamplingTest(ref, min, max, nullModel, x.asInstanceOf[Encode.VT])
  }

  def getStatistic(nm: NullModel, x: Encode.Coding): Double = {
    //println(s"scores: ${st.score.toArray.mkString(",")}")
    //println(s"variances: ${diag(st.variance).toArray.mkString(",")}")
    val m = x.asInstanceOf[Encode.VT].coding
    val ts = m.map{sv =>
      val st = ScoreTest(nm, sv)
      st.score(0)/st.variance(0, 0).sqrt
    }
    //val ts = st.score :/ diag(st.variance).map(x => x.sqrt)
    max(ts)
  }

  @SerialVersionUID(7727880101L)
  final case class AnalyticTest(nullModel: NullModel,
                                x: Encode.VT)
    extends VT with AssocMethod.AnalyticTest
  {
    val scoreTest = {
      //val cnt = x.coding.activeSize
      //println(s"geno: ${geno.coding.rows} x ${geno.coding.cols}, " +
      //  s"active: $cnt sample: ${geno.coding.toDense(0,::).t.toArray.mkString(",")}")
      val m = DenseVector.horzcat(x.coding:_*)
      ScoreTest(nullModel, m)
    }
    val statistic = getStatistic(nullModel, x)
    private val ss = diag(scoreTest.variance).map(_.sqrt)
    private val ts = scoreTest.score :/ ss
    val pValue = {
      val maxT = max(ts)
      val cutoff = maxT * DenseVector.ones[Double](ts.length)
      val pTry = Try(MultivariateNormal.Centered(scoreTest.variance).cdf(cutoff).pvalue)
      pTry match {
        case Success(p) => Some(1.0 - p)
        case _ => None
      }
    }
    def result: AssocMethod.VTAnalytic = {
      val info = s"Ts=${ts.toArray.mkString(",")}"
      AssocMethod.VTAnalytic(x.vars, x.size, statistic, pValue, info)
    }
  }

  @SerialVersionUID(7727880201L)
  final case class ResamplingTest(refStatistic: Double,
                                  min: Int,
                                  max: Int,
                                  nullModel: NullModel,
                                  x: Encode.VT)
    extends VT with AssocMethod.ResamplingTest
  {
    def pCount = {
      Resampling.Simple(refStatistic, min, max, nullModel, x, getStatistic).pCount
    }
    def result: AssocMethod.VTResampling =
      AssocMethod.VTResampling(x.vars, x.size, refStatistic, pCount)
  }

}

{'fraction_non_alphanumeric': 0.12715008431703204, 'fraction_numerical': 0.012141652613827993, 'mean_word_length': 2.601092896174863, 'pattern_counts': {'":': 0, '<': 3, '<?xml version=': 0, '>': 11, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 10, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'dir': 'git-commits-cleaned', 'id': '4549360', 'n_tokens_mistral': 2119, 'n_tokens_neox': 2015, 'n_words': 398}