diff --git a/calc.go b/calc.go
index e3db8826de41472bb35fdc4f804cb1a7ad93231d..ada84963fe35bfc9a0efdc45024295586274341b 100644
--- a/calc.go
+++ b/calc.go
@@ -333,6 +333,8 @@ type formulaFuncs struct {
 //    BESSELJ
 //    BESSELK
 //    BESSELY
+//    BETAINV
+//    BETA.INV
 //    BIN2DEC
 //    BIN2HEX
 //    BIN2OCT
@@ -415,6 +417,7 @@ type formulaFuncs struct {
 //    FALSE
 //    FIND
 //    FINDB
+//    F.INV.RT
 //    FINV
 //    FISHER
 //    FISHERINV
@@ -5187,966 +5190,1028 @@ func (fn *formulaFuncs) AVERAGEIF(argsList *list.List) formulaArg {
 	return newNumberFormulaArg(sum / count)
 }
 
-// incompleteGamma is an implementation of the incomplete gamma function.
-func incompleteGamma(a, x float64) float64 {
-	max := 32
-	summer := 0.0
-	for n := 0; n <= max; n++ {
-		divisor := a
-		for i := 1; i <= n; i++ {
-			divisor *= (a + float64(i))
-		}
-		summer += math.Pow(x, float64(n)) / divisor
+// d1mach returns double precision real machine constants.
+func d1mach(i int) float64 {
+	arr := []float64{
+		2.2250738585072014e-308,
+		1.7976931348623158e+308,
+		1.1102230246251565e-16,
+		2.2204460492503131e-16,
+		0.301029995663981195,
 	}
-	return math.Pow(x, a) * math.Exp(0-x) * summer
+	if i > len(arr) {
+		return 0
+	}
+	return arr[i-1]
 }
 
-// CHIDIST function calculates the right-tailed probability of the chi-square
-// distribution. The syntax of the function is:
-//
-//    CHIDIST(x,degrees_freedom)
-//
-func (fn *formulaFuncs) CHIDIST(argsList *list.List) formulaArg {
-	if argsList.Len() != 2 {
-		return newErrorFormulaArg(formulaErrorVALUE, "CHIDIST requires 2 numeric arguments")
-	}
-	x := argsList.Front().Value.(formulaArg).ToNumber()
-	if x.Type != ArgNumber {
-		return x
+// chebyshevInit determines the number of terms for the double precision
+// orthogonal series "dos" needed to insure the error is no larger
+// than "eta". Ordinarily eta will be chosen to be one-tenth machine
+// precision.
+func chebyshevInit(nos int, eta float64, dos []float64) int {
+	i, e := 0, 0.0
+	if nos < 1 {
+		return 0
 	}
-	degress := argsList.Back().Value.(formulaArg).ToNumber()
-	if degress.Type != ArgNumber {
-		return degress
+	for ii := 1; ii <= nos; ii++ {
+		i = nos - ii
+		e += math.Abs(dos[i])
+		if e > eta {
+			return i
+		}
 	}
-	return newNumberFormulaArg(1 - (incompleteGamma(degress.Number/2, x.Number/2) / math.Gamma(degress.Number/2)))
+	return i
 }
 
-// confidence is an implementation of the formula functions CONFIDENCE and
-// CONFIDENCE.NORM.
-func (fn *formulaFuncs) confidence(name string, argsList *list.List) formulaArg {
-	if argsList.Len() != 3 {
-		return newErrorFormulaArg(formulaErrorVALUE, fmt.Sprintf("%s requires 3 numeric arguments", name))
+// chebyshevEval evaluates the n-term Chebyshev series "a" at "x".
+func chebyshevEval(n int, x float64, a []float64) float64 {
+	if n < 1 || n > 1000 || x < -1.1 || x > 1.1 {
+		return math.NaN()
 	}
-	alpha := argsList.Front().Value.(formulaArg).ToNumber()
-	if alpha.Type != ArgNumber {
-		return alpha
+	twox, b0, b1, b2 := x*2, 0.0, 0.0, 0.0
+	for i := 1; i <= n; i++ {
+		b2 = b1
+		b1 = b0
+		b0 = twox*b1 - b2 + a[n-i]
 	}
-	if alpha.Number <= 0 || alpha.Number >= 1 {
-		return newErrorFormulaArg(formulaErrorNUM, formulaErrorNUM)
+	return (b0 - b2) * 0.5
+}
+
+// lgammacor is an implementation for the log(gamma) correction.
+func lgammacor(x float64) float64 {
+	algmcs := []float64{
+		0.1666389480451863247205729650822, -0.1384948176067563840732986059135e-4,
+		0.9810825646924729426157171547487e-8, -0.1809129475572494194263306266719e-10,
+		0.6221098041892605227126015543416e-13, -0.3399615005417721944303330599666e-15,
+		0.2683181998482698748957538846666e-17, -0.2868042435334643284144622399999e-19,
+		0.3962837061046434803679306666666e-21, -0.6831888753985766870111999999999e-23,
+		0.1429227355942498147573333333333e-24, -0.3547598158101070547199999999999e-26,
+		0.1025680058010470912000000000000e-27, -0.3401102254316748799999999999999e-29,
+		0.1276642195630062933333333333333e-30,
 	}
-	stdDev := argsList.Front().Next().Value.(formulaArg).ToNumber()
-	if stdDev.Type != ArgNumber {
-		return stdDev
+	nalgm := chebyshevInit(15, d1mach(3), algmcs)
+	xbig := 1.0 / math.Sqrt(d1mach(3))
+	xmax := math.Exp(math.Min(math.Log(d1mach(2)/12.0), -math.Log(12.0*d1mach(1))))
+	if x < 10.0 {
+		return math.NaN()
+	} else if x >= xmax {
+		return 4.930380657631324e-32
+	} else if x < xbig {
+		tmp := 10.0 / x
+		return chebyshevEval(nalgm, tmp*tmp*2.0-1.0, algmcs) / x
 	}
-	if stdDev.Number <= 0 {
-		return newErrorFormulaArg(formulaErrorNUM, formulaErrorNUM)
+	return 1.0 / (x * 12.0)
+}
+
+// logrelerr compute the relative error logarithm.
+func logrelerr(x float64) float64 {
+	alnrcs := []float64{
+		0.10378693562743769800686267719098e+1, -0.13364301504908918098766041553133,
+		0.19408249135520563357926199374750e-1, -0.30107551127535777690376537776592e-2,
+		0.48694614797154850090456366509137e-3, -0.81054881893175356066809943008622e-4,
+		0.13778847799559524782938251496059e-4, -0.23802210894358970251369992914935e-5,
+		0.41640416213865183476391859901989e-6, -0.73595828378075994984266837031998e-7,
+		0.13117611876241674949152294345011e-7, -0.23546709317742425136696092330175e-8,
+		0.42522773276034997775638052962567e-9, -0.77190894134840796826108107493300e-10,
+		0.14075746481359069909215356472191e-10, -0.25769072058024680627537078627584e-11,
+		0.47342406666294421849154395005938e-12, -0.87249012674742641745301263292675e-13,
+		0.16124614902740551465739833119115e-13, -0.29875652015665773006710792416815e-14,
+		0.55480701209082887983041321697279e-15, -0.10324619158271569595141333961932e-15,
+		0.19250239203049851177878503244868e-16, -0.35955073465265150011189707844266e-17,
+		0.67264542537876857892194574226773e-18, -0.12602624168735219252082425637546e-18,
+		0.23644884408606210044916158955519e-19, -0.44419377050807936898878389179733e-20,
+		0.83546594464034259016241293994666e-21, -0.15731559416479562574899253521066e-21,
+		0.29653128740247422686154369706666e-22, -0.55949583481815947292156013226666e-23,
+		0.10566354268835681048187284138666e-23, -0.19972483680670204548314999466666e-24,
+		0.37782977818839361421049855999999e-25, -0.71531586889081740345038165333333e-26,
+		0.13552488463674213646502024533333e-26, -0.25694673048487567430079829333333e-27,
+		0.48747756066216949076459519999999e-28, -0.92542112530849715321132373333333e-29,
+		0.17578597841760239233269760000000e-29, -0.33410026677731010351377066666666e-30,
+		0.63533936180236187354180266666666e-31,
 	}
-	size := argsList.Back().Value.(formulaArg).ToNumber()
-	if size.Type != ArgNumber {
-		return size
+	nlnrel := chebyshevInit(43, 0.1*d1mach(3), alnrcs)
+	if x <= -1 {
+		return math.NaN()
 	}
-	if size.Number < 1 {
-		return newErrorFormulaArg(formulaErrorNUM, formulaErrorNUM)
+	if math.Abs(x) <= 0.375 {
+		return x * (1.0 - x*chebyshevEval(nlnrel, x/0.375, alnrcs))
 	}
-	args := list.New()
-	args.Init()
-	args.PushBack(newNumberFormulaArg(alpha.Number / 2))
-	args.PushBack(newNumberFormulaArg(0))
-	args.PushBack(newNumberFormulaArg(1))
-	return newNumberFormulaArg(-fn.NORMINV(args).Number * (stdDev.Number / math.Sqrt(size.Number)))
-}
-
-// CONFIDENCE function uses a Normal Distribution to calculate a confidence
-// value that can be used to construct the Confidence Interval for a
-// population mean, for a supplied probablity and sample size. It is assumed
-// that the standard deviation of the population is known. The syntax of the
-// function is:
-//
-//    CONFIDENCE(alpha,standard_dev,size)
-//
-func (fn *formulaFuncs) CONFIDENCE(argsList *list.List) formulaArg {
-	return fn.confidence("CONFIDENCE", argsList)
-}
-
-// CONFIDENCEdotNORM function uses a Normal Distribution to calculate a
-// confidence value that can be used to construct the confidence interval for
-// a population mean, for a supplied probablity and sample size. It is
-// assumed that the standard deviation of the population is known. The syntax
-// of the Confidence.Norm function is:
-//
-//    CONFIDENCE.NORM(alpha,standard_dev,size)
-//
-func (fn *formulaFuncs) CONFIDENCEdotNORM(argsList *list.List) formulaArg {
-	return fn.confidence("CONFIDENCE.NORM", argsList)
+	return math.Log(x + 1.0)
 }
 
-// COVAR function calculates the covariance of two supplied sets of values. The
-// syntax of the function is:
-//
-//    COVAR(array1,array2)
-//
-func (fn *formulaFuncs) COVAR(argsList *list.List) formulaArg {
-	if argsList.Len() != 2 {
-		return newErrorFormulaArg(formulaErrorVALUE, "COVAR requires 2 arguments")
+// logBeta is an implementation for the log of the beta distribution
+// function.
+func logBeta(a, b float64) float64 {
+	corr, p, q := 0.0, a, a
+	if b < p {
+		p = b
 	}
-	array1 := argsList.Front().Value.(formulaArg)
-	array2 := argsList.Back().Value.(formulaArg)
-	left, right := array1.ToList(), array2.ToList()
-	n := len(left)
-	if n != len(right) {
-		return newErrorFormulaArg(formulaErrorNA, formulaErrorNA)
+	if b > q {
+		q = b
 	}
-	l1, l2 := list.New(), list.New()
-	l1.PushBack(array1)
-	l2.PushBack(array2)
-	result, skip := 0.0, 0
-	mean1, mean2 := fn.AVERAGE(l1), fn.AVERAGE(l2)
-	for i := 0; i < n; i++ {
-		arg1 := left[i].ToNumber()
-		arg2 := right[i].ToNumber()
-		if arg1.Type == ArgError || arg2.Type == ArgError {
-			skip++
-			continue
-		}
-		result += (arg1.Number - mean1.Number) * (arg2.Number - mean2.Number)
+	if p < 0 {
+		return math.NaN()
 	}
-	return newNumberFormulaArg(result / float64(n-skip))
-}
-
-// COVARIANCEdotP function calculates the population covariance of two supplied
-// sets of values. The syntax of the function is:
-//
-//    COVARIANCE.P(array1,array2)
-//
-func (fn *formulaFuncs) COVARIANCEdotP(argsList *list.List) formulaArg {
-	if argsList.Len() != 2 {
-		return newErrorFormulaArg(formulaErrorVALUE, "COVARIANCE.P requires 2 arguments")
+	if p == 0 {
+		return math.MaxFloat64
 	}
-	return fn.COVAR(argsList)
-}
-
-// calcStringCountSum is part of the implementation countSum.
-func calcStringCountSum(countText bool, count, sum float64, num, arg formulaArg) (float64, float64) {
-	if countText && num.Type == ArgError && arg.String != "" {
-		count++
+	if p >= 10.0 {
+		corr = lgammacor(p) + lgammacor(q) - lgammacor(p+q)
+		return math.Log(q)*-0.5 + 0.918938533204672741780329736406 + corr + (p-0.5)*math.Log(p/(p+q)) + q*logrelerr(-p/(p+q))
 	}
-	if num.Type == ArgNumber {
-		sum += num.Number
-		count++
+	if q >= 10 {
+		corr = lgammacor(q) - lgammacor(p+q)
+		val, _ := math.Lgamma(p)
+		return val + corr + p - p*math.Log(p+q) + (q-0.5)*logrelerr(-p/(p+q))
 	}
-	return count, sum
+	return math.Log(math.Gamma(p) * (math.Gamma(q) / math.Gamma(p+q)))
 }
 
-// countSum get count and sum for a formula arguments array.
-func (fn *formulaFuncs) countSum(countText bool, args []formulaArg) (count, sum float64) {
-	for _, arg := range args {
-		switch arg.Type {
-		case ArgNumber:
-			if countText || !arg.Boolean {
-				sum += arg.Number
-				count++
-			}
-		case ArgString:
-			if !countText && (arg.Value() == "TRUE" || arg.Value() == "FALSE") {
-				continue
-			} else if countText && (arg.Value() == "TRUE" || arg.Value() == "FALSE") {
-				num := arg.ToBool()
-				if num.Type == ArgNumber {
-					count++
-					sum += num.Number
-					continue
-				}
+// pbetaRaw is a part of pbeta for the beta distribution.
+func pbetaRaw(alnsml, ans, eps, p, pin, q, sml, x, y float64) float64 {
+	if q > 1.0 {
+		xb := p*math.Log(y) + q*math.Log(1.0-y) - logBeta(p, q) - math.Log(q)
+		ib := int(math.Max(xb/alnsml, 0.0))
+		term := math.Exp(xb - float64(ib)*alnsml)
+		c := 1.0 / (1.0 - y)
+		p1 := q * c / (p + q - 1.0)
+		finsum := 0.0
+		n := int(q)
+		if q == float64(n) {
+			n = n - 1
+		}
+		for i := 1; i <= n; i++ {
+			if p1 <= 1 && term/eps <= finsum {
+				break
+			}
+			xi := float64(i)
+			term = (q - xi + 1.0) * c * term / (p + q - xi)
+			if term > 1.0 {
+				ib = ib - 1
+				term = term * sml
+			}
+			if ib == 0 {
+				finsum = finsum + term
 			}
-			num := arg.ToNumber()
-			count, sum = calcStringCountSum(countText, count, sum, num, arg)
-		case ArgList, ArgMatrix:
-			cnt, summary := fn.countSum(countText, arg.ToList())
-			sum += summary
-			count += cnt
 		}
+		ans = ans + finsum
 	}
-	return
+	if y != x || p != pin {
+		ans = 1.0 - ans
+	}
+	ans = math.Max(math.Min(ans, 1.0), 0.0)
+	return ans
 }
 
-// CORREL function calculates the Pearson Product-Moment Correlation
-// Coefficient for two sets of values. The syntax of the function is:
-//
-//    CORREL(array1,array2)
-//
-func (fn *formulaFuncs) CORREL(argsList *list.List) formulaArg {
-	if argsList.Len() != 2 {
-		return newErrorFormulaArg(formulaErrorVALUE, "CORREL requires 2 arguments")
-	}
-	array1 := argsList.Front().Value.(formulaArg)
-	array2 := argsList.Back().Value.(formulaArg)
-	left, right := array1.ToList(), array2.ToList()
-	n := len(left)
-	if n != len(right) {
-		return newErrorFormulaArg(formulaErrorNA, formulaErrorNA)
+// pbeta returns distribution function of the beta distribution.
+func pbeta(x, pin, qin float64) (ans float64) {
+	eps := d1mach(3)
+	alneps := math.Log(eps)
+	sml := d1mach(1)
+	alnsml := math.Log(sml)
+	y := x
+	p := pin
+	q := qin
+	if p/(p+q) < x {
+		y = 1.0 - y
+		p = qin
+		q = pin
 	}
-	l1, l2, l3 := list.New(), list.New(), list.New()
-	for i := 0; i < n; i++ {
-		if lhs, rhs := left[i].ToNumber(), right[i].ToNumber(); lhs.Number != 0 && rhs.Number != 0 {
-			l1.PushBack(lhs)
-			l2.PushBack(rhs)
+	if (p+q)*y/(p+1.0) < eps {
+		xb := p*math.Log(math.Max(y, sml)) - math.Log(p) - logBeta(p, q)
+		if xb > alnsml && y != 0.0 {
+			ans = math.Exp(xb)
 		}
-	}
-	stdev1, stdev2 := fn.STDEV(l1), fn.STDEV(l2)
-	if stdev1.Number == 0 || stdev2.Number == 0 {
-		return newErrorFormulaArg(formulaErrorDIV, formulaErrorDIV)
-	}
-	mean1, mean2, skip := fn.AVERAGE(l1), fn.AVERAGE(l2), 0
-	for i := 0; i < n; i++ {
-		lhs, rhs := left[i].ToNumber(), right[i].ToNumber()
-		if lhs.Number == 0 || rhs.Number == 0 {
-			skip++
-			continue
+		if y != x || p != pin {
+			ans = 1.0 - ans
 		}
-		l3.PushBack(newNumberFormulaArg((lhs.Number - mean1.Number) * (rhs.Number - mean2.Number)))
-	}
-	return newNumberFormulaArg(fn.SUM(l3).Number / float64(n-skip-1) / stdev1.Number / stdev2.Number)
-}
-
-// COUNT function returns the count of numeric values in a supplied set of
-// cells or values. This count includes both numbers and dates. The syntax of
-// the function is:
-//
-//    COUNT(value1,[value2],...)
-//
-func (fn *formulaFuncs) COUNT(argsList *list.List) formulaArg {
-	var count int
-	for token := argsList.Front(); token != nil; token = token.Next() {
-		arg := token.Value.(formulaArg)
-		switch arg.Type {
-		case ArgString:
-			if arg.ToNumber().Type != ArgError {
-				count++
-			}
-		case ArgNumber:
-			count++
-		case ArgMatrix:
-			for _, row := range arg.Matrix {
-				for _, value := range row {
-					if value.ToNumber().Type != ArgError {
-						count++
-					}
+	} else {
+		ps := q - math.Floor(q)
+		if ps == 0.0 {
+			ps = 1.0
+		}
+		xb := p*math.Log(y) - logBeta(ps, p) - math.Log(p)
+		if xb >= alnsml {
+			ans = math.Exp(xb)
+			term := ans * p
+			if ps != 1.0 {
+				n := int(math.Max(alneps/math.Log(y), 4.0))
+				for i := 1; i <= n; i++ {
+					xi := float64(i)
+					term = term * (xi - ps) * y / xi
+					ans = ans + term/(p+xi)
 				}
 			}
 		}
+		ans = pbetaRaw(alnsml, ans, eps, p, pin, q, sml, x, y)
 	}
-	return newNumberFormulaArg(float64(count))
+	return ans
 }
 
-// COUNTA function returns the number of non-blanks within a supplied set of
-// cells or values. The syntax of the function is:
-//
-//    COUNTA(value1,[value2],...)
-//
-func (fn *formulaFuncs) COUNTA(argsList *list.List) formulaArg {
-	var count int
-	for token := argsList.Front(); token != nil; token = token.Next() {
-		arg := token.Value.(formulaArg)
-		switch arg.Type {
-		case ArgString:
-			if arg.String != "" {
-				count++
-			}
-		case ArgNumber:
-			count++
-		case ArgMatrix:
-			for _, row := range arg.ToList() {
-				switch row.Type {
-				case ArgString:
-					if row.String != "" {
-						count++
+// betainvProbIterator is a part of betainv for the inverse of the beta
+// function.
+func betainvProbIterator(alpha1, alpha3, beta1, beta2, beta3, logbeta, lower, maxCumulative, prob1, prob2, upper float64, needSwap bool) float64 {
+	var i, j, prev, prop4 float64
+	j = 1
+	for prob := 0; prob < 1000; prob++ {
+		prop3 := pbeta(beta3, alpha1, beta1)
+		prop3 = (prop3 - prob1) * math.Exp(logbeta+prob2*math.Log(beta3)+beta2*math.Log(1.0-beta3))
+		if prop3*prop4 <= 0 {
+			prev = math.Max(math.Abs(j), maxCumulative)
+		}
+		h := 1.0
+		for iteratorCount := 0; iteratorCount < 1000; iteratorCount++ {
+			j = h * prop3
+			if math.Abs(j) < prev {
+				i = beta3 - j
+				if i >= 0 && i <= 1.0 {
+					if prev <= alpha3 {
+						return beta3
+					}
+					if math.Abs(prop3) <= alpha3 {
+						return beta3
+					}
+					if i != 0 && i != 1.0 {
+						break
 					}
-				case ArgNumber:
-					count++
 				}
 			}
+			h /= 3.0
+		}
+		if i == beta3 {
+			return beta3
 		}
+		beta3, prop4 = i, prop3
 	}
-	return newNumberFormulaArg(float64(count))
+	return beta3
 }
 
-// COUNTBLANK function returns the number of blank cells in a supplied range.
-// The syntax of the function is:
-//
-//    COUNTBLANK(range)
-//
-func (fn *formulaFuncs) COUNTBLANK(argsList *list.List) formulaArg {
-	if argsList.Len() != 1 {
-		return newErrorFormulaArg(formulaErrorVALUE, "COUNTBLANK requires 1 argument")
+// calcBetainv is an implementation for the quantile of the beta
+// distribution.
+func calcBetainv(probability, alpha, beta, lower, upper float64) float64 {
+	minCumulative, maxCumulative := 1.0e-300, 3.0e-308
+	lowerBound, upperBound := maxCumulative, 1.0-2.22e-16
+	needSwap := false
+	var alpha1, alpha2, beta1, beta2, beta3, prob1, x, y float64
+	if probability <= 0.5 {
+		prob1, alpha1, beta1 = probability, alpha, beta
+	} else {
+		prob1, alpha1, beta1, needSwap = 1.0-probability, beta, alpha, true
 	}
-	var count float64
-	for _, cell := range argsList.Front().Value.(formulaArg).ToList() {
-		if cell.Value() == "" {
-			count++
+	logbeta := logBeta(alpha, beta)
+	prob2 := math.Sqrt(-math.Log(prob1 * prob1))
+	prob3 := prob2 - (prob2*0.27061+2.3075)/(prob2*(prob2*0.04481+0.99229)+1)
+	if alpha1 > 1 && beta1 > 1 {
+		alpha2, beta2, prob2 = 1/(alpha1+alpha1-1), 1/(beta1+beta1-1), (prob3*prob3-3)/6
+		x = 2 / (alpha2 + beta2)
+		y = prob3*math.Sqrt(x+prob2)/x - (beta2-alpha2)*(prob2+5/6.0-2/(x*3))
+		beta3 = alpha1 / (alpha1 + beta1*math.Exp(y+y))
+	} else {
+		beta2, prob2 = 1/(beta1*9), beta1+beta1
+		beta2 = prob2 * math.Pow(1-beta2+prob3*math.Sqrt(beta2), 3)
+		if beta2 <= 0 {
+			beta3 = 1 - math.Exp((math.Log((1-prob1)*beta1)+logbeta)/beta1)
+		} else {
+			beta2 = (prob2 + alpha1*4 - 2) / beta2
+			if beta2 <= 1 {
+				beta3 = math.Exp((logbeta + math.Log(alpha1*prob1)) / alpha1)
+			} else {
+				beta3 = 1 - 2/(beta2+1)
+			}
 		}
 	}
-	return newNumberFormulaArg(count)
+	beta2, prob2 = 1-beta1, 1-alpha1
+	if beta3 < lowerBound {
+		beta3 = lowerBound
+	} else if beta3 > upperBound {
+		beta3 = upperBound
+	}
+	alpha3 := math.Max(minCumulative, math.Pow(10.0, -13.0-2.5/(alpha1*alpha1)-0.5/(prob1*prob1)))
+	beta3 = betainvProbIterator(alpha1, alpha3, beta1, beta2, beta3, logbeta, lower, maxCumulative, prob1, prob2, upper, needSwap)
+	if needSwap {
+		beta3 = 1.0 - beta3
+	}
+	return (upper-lower)*beta3 + lower
 }
 
-// COUNTIF function returns the number of cells within a supplied range, that
-// satisfy a given criteria. The syntax of the function is:
-//
-//    COUNTIF(range,criteria)
-//
-func (fn *formulaFuncs) COUNTIF(argsList *list.List) formulaArg {
-	if argsList.Len() != 2 {
-		return newErrorFormulaArg(formulaErrorVALUE, "COUNTIF requires 2 arguments")
+// betainv is an implementation of the formula functions BETAINV and
+// BETA.INV.
+func (fn *formulaFuncs) betainv(name string, argsList *list.List) formulaArg {
+	if argsList.Len() < 3 {
+		return newErrorFormulaArg(formulaErrorVALUE, fmt.Sprintf("%s requires at least 3 arguments", name))
 	}
-	var (
-		criteria = formulaCriteriaParser(argsList.Front().Next().Value.(formulaArg).String)
-		count    float64
-	)
-	for _, cell := range argsList.Front().Value.(formulaArg).ToList() {
-		if ok, _ := formulaCriteriaEval(cell.Value(), criteria); ok {
-			count++
-		}
+	if argsList.Len() > 5 {
+		return newErrorFormulaArg(formulaErrorVALUE, fmt.Sprintf("%s requires at most 5 arguments", name))
 	}
-	return newNumberFormulaArg(count)
-}
-
-// formulaIfsMatch function returns cells reference array which match criteria.
-func formulaIfsMatch(args []formulaArg) (cellRefs []cellRef) {
-	for i := 0; i < len(args)-1; i += 2 {
-		match := []cellRef{}
-		matrix, criteria := args[i].Matrix, formulaCriteriaParser(args[i+1].Value())
-		if i == 0 {
-			for rowIdx, row := range matrix {
-				for colIdx, col := range row {
-					if ok, _ := formulaCriteriaEval(col.Value(), criteria); ok {
-						match = append(match, cellRef{Col: colIdx, Row: rowIdx})
-					}
-				}
-			}
-		} else {
-			for _, ref := range cellRefs {
-				value := matrix[ref.Row][ref.Col]
-				if ok, _ := formulaCriteriaEval(value.Value(), criteria); ok {
-					match = append(match, ref)
-				}
-			}
-		}
-		if len(match) == 0 {
-			return
-		}
-		cellRefs = match[:]
+	probability := argsList.Front().Value.(formulaArg).ToNumber()
+	if probability.Type != ArgNumber {
+		return probability
 	}
-	return
-}
-
-// COUNTIFS function returns the number of rows within a table, that satisfy a
-// set of given criteria. The syntax of the function is:
-//
-//    COUNTIFS(criteria_range1,criteria1,[criteria_range2,criteria2],...)
-//
-func (fn *formulaFuncs) COUNTIFS(argsList *list.List) formulaArg {
-	if argsList.Len() < 2 {
-		return newErrorFormulaArg(formulaErrorVALUE, "COUNTIFS requires at least 2 arguments")
+	if probability.Number <= 0 || probability.Number >= 1 {
+		return newErrorFormulaArg(formulaErrorNUM, formulaErrorNUM)
 	}
-	if argsList.Len()%2 != 0 {
-		return newErrorFormulaArg(formulaErrorNA, formulaErrorNA)
+	alpha := argsList.Front().Next().Value.(formulaArg).ToNumber()
+	if alpha.Type != ArgNumber {
+		return alpha
 	}
-	args := []formulaArg{}
-	for arg := argsList.Front(); arg != nil; arg = arg.Next() {
-		args = append(args, arg.Value.(formulaArg))
+	beta := argsList.Front().Next().Next().Value.(formulaArg).ToNumber()
+	if beta.Type != ArgNumber {
+		return beta
 	}
-	return newNumberFormulaArg(float64(len(formulaIfsMatch(args))))
-}
-
-// DEVSQ function calculates the sum of the squared deviations from the sample
-// mean. The syntax of the function is:
-//
-//    DEVSQ(number1,[number2],...)
-//
-func (fn *formulaFuncs) DEVSQ(argsList *list.List) formulaArg {
-	if argsList.Len() < 1 {
-		return newErrorFormulaArg(formulaErrorVALUE, "DEVSQ requires at least 1 numeric argument")
+	if alpha.Number <= 0 || beta.Number <= 0 {
+		return newErrorFormulaArg(formulaErrorNUM, formulaErrorNUM)
 	}
-	avg, count, result := fn.AVERAGE(argsList), -1, 0.0
-	for arg := argsList.Front(); arg != nil; arg = arg.Next() {
-		for _, number := range arg.Value.(formulaArg).ToList() {
-			num := number.ToNumber()
-			if num.Type != ArgNumber {
-				continue
-			}
-			count++
-			if count == 0 {
-				result = math.Pow(num.Number-avg.Number, 2)
-				continue
-			}
-			result += math.Pow(num.Number-avg.Number, 2)
+	a, b := newNumberFormulaArg(0), newNumberFormulaArg(1)
+	if argsList.Len() > 3 {
+		if a = argsList.Front().Next().Next().Next().Value.(formulaArg).ToNumber(); a.Type != ArgNumber {
+			return a
 		}
 	}
-	if count == -1 {
-		return newErrorFormulaArg(formulaErrorNA, formulaErrorNA)
+	if argsList.Len() == 5 {
+		if b = argsList.Back().Value.(formulaArg).ToNumber(); b.Type != ArgNumber {
+			return b
+		}
 	}
-	return newNumberFormulaArg(result)
+	if a.Number == b.Number {
+		return newErrorFormulaArg(formulaErrorNUM, formulaErrorNUM)
+	}
+	return newNumberFormulaArg(calcBetainv(probability.Number, alpha.Number, beta.Number, a.Number, b.Number))
 }
 
-// FISHER function calculates the Fisher Transformation for a supplied value.
-// The syntax of the function is:
+// BETAINV function uses an iterative procedure to calculate the inverse of
+// the cumulative beta probability density function for a supplied
+// probability. The syntax of the function is:
 //
-//    FISHER(x)
+//    BETAINV(probability,alpha,beta,[A],[B])
 //
-func (fn *formulaFuncs) FISHER(argsList *list.List) formulaArg {
-	if argsList.Len() != 1 {
-		return newErrorFormulaArg(formulaErrorVALUE, "FISHER requires 1 numeric argument")
-	}
-	token := argsList.Front().Value.(formulaArg)
-	switch token.Type {
-	case ArgString:
-		arg := token.ToNumber()
-		if arg.Type == ArgNumber {
-			if arg.Number <= -1 || arg.Number >= 1 {
-				return newErrorFormulaArg(formulaErrorNA, formulaErrorNA)
-			}
-			return newNumberFormulaArg(0.5 * math.Log((1+arg.Number)/(1-arg.Number)))
-		}
-	case ArgNumber:
-		if token.Number <= -1 || token.Number >= 1 {
-			return newErrorFormulaArg(formulaErrorNA, formulaErrorNA)
-		}
-		return newNumberFormulaArg(0.5 * math.Log((1+token.Number)/(1-token.Number)))
-	}
-	return newErrorFormulaArg(formulaErrorVALUE, "FISHER requires 1 numeric argument")
+func (fn *formulaFuncs) BETAINV(argsList *list.List) formulaArg {
+	return fn.betainv("BETAINV", argsList)
 }
 
-// FISHERINV function calculates the inverse of the Fisher Transformation and
-// returns a value between -1 and +1. The syntax of the function is:
+// BETAdotINV function uses an iterative procedure to calculate the inverse of
+// the cumulative beta probability density function for a supplied
+// probability. The syntax of the function is:
 //
-//    FISHERINV(y)
+//    BETA.INV(probability,alpha,beta,[A],[B])
 //
-func (fn *formulaFuncs) FISHERINV(argsList *list.List) formulaArg {
-	if argsList.Len() != 1 {
-		return newErrorFormulaArg(formulaErrorVALUE, "FISHERINV requires 1 numeric argument")
-	}
-	token := argsList.Front().Value.(formulaArg)
-	switch token.Type {
-	case ArgString:
-		arg := token.ToNumber()
-		if arg.Type == ArgNumber {
-			return newNumberFormulaArg((math.Exp(2*arg.Number) - 1) / (math.Exp(2*arg.Number) + 1))
-		}
-	case ArgNumber:
-		return newNumberFormulaArg((math.Exp(2*token.Number) - 1) / (math.Exp(2*token.Number) + 1))
-	}
-	return newErrorFormulaArg(formulaErrorVALUE, "FISHERINV requires 1 numeric argument")
+func (fn *formulaFuncs) BETAdotINV(argsList *list.List) formulaArg {
+	return fn.betainv("BETA.INV", argsList)
 }
 
-// GAMMA function returns the value of the Gamma Function, Γ(n), for a
-// specified number, n. The syntax of the function is:
-//
-//    GAMMA(number)
-//
-func (fn *formulaFuncs) GAMMA(argsList *list.List) formulaArg {
-	if argsList.Len() != 1 {
-		return newErrorFormulaArg(formulaErrorVALUE, "GAMMA requires 1 numeric argument")
-	}
-	token := argsList.Front().Value.(formulaArg)
-	switch token.Type {
-	case ArgString:
-		arg := token.ToNumber()
-		if arg.Type == ArgNumber {
-			if arg.Number <= 0 {
-				return newErrorFormulaArg(formulaErrorNA, formulaErrorNA)
-			}
-			return newNumberFormulaArg(math.Gamma(arg.Number))
-		}
-	case ArgNumber:
-		if token.Number <= 0 {
-			return newErrorFormulaArg(formulaErrorNA, formulaErrorNA)
+// incompleteGamma is an implementation of the incomplete gamma function.
+func incompleteGamma(a, x float64) float64 {
+	max := 32
+	summer := 0.0
+	for n := 0; n <= max; n++ {
+		divisor := a
+		for i := 1; i <= n; i++ {
+			divisor *= (a + float64(i))
 		}
-		return newNumberFormulaArg(math.Gamma(token.Number))
+		summer += math.Pow(x, float64(n)) / divisor
 	}
-	return newErrorFormulaArg(formulaErrorVALUE, "GAMMA requires 1 numeric argument")
+	return math.Pow(x, a) * math.Exp(0-x) * summer
 }
 
-// GAMMALN function returns the natural logarithm of the Gamma Function, Γ
-// (n). The syntax of the function is:
+// CHIDIST function calculates the right-tailed probability of the chi-square
+// distribution. The syntax of the function is:
 //
-//    GAMMALN(x)
+//    CHIDIST(x,degrees_freedom)
 //
-func (fn *formulaFuncs) GAMMALN(argsList *list.List) formulaArg {
-	if argsList.Len() != 1 {
-		return newErrorFormulaArg(formulaErrorVALUE, "GAMMALN requires 1 numeric argument")
+func (fn *formulaFuncs) CHIDIST(argsList *list.List) formulaArg {
+	if argsList.Len() != 2 {
+		return newErrorFormulaArg(formulaErrorVALUE, "CHIDIST requires 2 numeric arguments")
 	}
-	token := argsList.Front().Value.(formulaArg)
-	switch token.Type {
-	case ArgString:
-		arg := token.ToNumber()
-		if arg.Type == ArgNumber {
-			if arg.Number <= 0 {
-				return newErrorFormulaArg(formulaErrorNA, formulaErrorNA)
-			}
-			return newNumberFormulaArg(math.Log(math.Gamma(arg.Number)))
-		}
-	case ArgNumber:
-		if token.Number <= 0 {
-			return newErrorFormulaArg(formulaErrorNA, formulaErrorNA)
-		}
-		return newNumberFormulaArg(math.Log(math.Gamma(token.Number)))
+	x := argsList.Front().Value.(formulaArg).ToNumber()
+	if x.Type != ArgNumber {
+		return x
 	}
-	return newErrorFormulaArg(formulaErrorVALUE, "GAMMALN requires 1 numeric argument")
+	degress := argsList.Back().Value.(formulaArg).ToNumber()
+	if degress.Type != ArgNumber {
+		return degress
+	}
+	return newNumberFormulaArg(1 - (incompleteGamma(degress.Number/2, x.Number/2) / math.Gamma(degress.Number/2)))
 }
 
-// GEOMEAN function calculates the geometric mean of a supplied set of values.
-// The syntax of the function is:
-//
-//    GEOMEAN(number1,[number2],...)
-//
-func (fn *formulaFuncs) GEOMEAN(argsList *list.List) formulaArg {
-	if argsList.Len() < 1 {
-		return newErrorFormulaArg(formulaErrorVALUE, "GEOMEAN requires at least 1 numeric argument")
+// confidence is an implementation of the formula functions CONFIDENCE and
+// CONFIDENCE.NORM.
+func (fn *formulaFuncs) confidence(name string, argsList *list.List) formulaArg {
+	if argsList.Len() != 3 {
+		return newErrorFormulaArg(formulaErrorVALUE, fmt.Sprintf("%s requires 3 numeric arguments", name))
 	}
-	product := fn.PRODUCT(argsList)
-	if product.Type != ArgNumber {
-		return product
+	alpha := argsList.Front().Value.(formulaArg).ToNumber()
+	if alpha.Type != ArgNumber {
+		return alpha
 	}
-	count := fn.COUNT(argsList)
-	min := fn.MIN(argsList)
-	if product.Number > 0 && min.Number > 0 {
-		return newNumberFormulaArg(math.Pow(product.Number, (1 / count.Number)))
+	if alpha.Number <= 0 || alpha.Number >= 1 {
+		return newErrorFormulaArg(formulaErrorNUM, formulaErrorNUM)
 	}
-	return newErrorFormulaArg(formulaErrorNUM, formulaErrorNUM)
+	stdDev := argsList.Front().Next().Value.(formulaArg).ToNumber()
+	if stdDev.Type != ArgNumber {
+		return stdDev
+	}
+	if stdDev.Number <= 0 {
+		return newErrorFormulaArg(formulaErrorNUM, formulaErrorNUM)
+	}
+	size := argsList.Back().Value.(formulaArg).ToNumber()
+	if size.Type != ArgNumber {
+		return size
+	}
+	if size.Number < 1 {
+		return newErrorFormulaArg(formulaErrorNUM, formulaErrorNUM)
+	}
+	args := list.New()
+	args.Init()
+	args.PushBack(newNumberFormulaArg(alpha.Number / 2))
+	args.PushBack(newNumberFormulaArg(0))
+	args.PushBack(newNumberFormulaArg(1))
+	return newNumberFormulaArg(-fn.NORMINV(args).Number * (stdDev.Number / math.Sqrt(size.Number)))
 }
 
-// HARMEAN function calculates the harmonic mean of a supplied set of values.
-// The syntax of the function is:
+// CONFIDENCE function uses a Normal Distribution to calculate a confidence
+// value that can be used to construct the Confidence Interval for a
+// population mean, for a supplied probablity and sample size. It is assumed
+// that the standard deviation of the population is known. The syntax of the
+// function is:
 //
-//    HARMEAN(number1,[number2],...)
+//    CONFIDENCE(alpha,standard_dev,size)
 //
-func (fn *formulaFuncs) HARMEAN(argsList *list.List) formulaArg {
-	if argsList.Len() < 1 {
-		return newErrorFormulaArg(formulaErrorVALUE, "HARMEAN requires at least 1 argument")
+func (fn *formulaFuncs) CONFIDENCE(argsList *list.List) formulaArg {
+	return fn.confidence("CONFIDENCE", argsList)
+}
+
+// CONFIDENCEdotNORM function uses a Normal Distribution to calculate a
+// confidence value that can be used to construct the confidence interval for
+// a population mean, for a supplied probablity and sample size. It is
+// assumed that the standard deviation of the population is known. The syntax
+// of the Confidence.Norm function is:
+//
+//    CONFIDENCE.NORM(alpha,standard_dev,size)
+//
+func (fn *formulaFuncs) CONFIDENCEdotNORM(argsList *list.List) formulaArg {
+	return fn.confidence("CONFIDENCE.NORM", argsList)
+}
+
+// COVAR function calculates the covariance of two supplied sets of values. The
+// syntax of the function is:
+//
+//    COVAR(array1,array2)
+//
+func (fn *formulaFuncs) COVAR(argsList *list.List) formulaArg {
+	if argsList.Len() != 2 {
+		return newErrorFormulaArg(formulaErrorVALUE, "COVAR requires 2 arguments")
 	}
-	if min := fn.MIN(argsList); min.Number < 0 {
+	array1 := argsList.Front().Value.(formulaArg)
+	array2 := argsList.Back().Value.(formulaArg)
+	left, right := array1.ToList(), array2.ToList()
+	n := len(left)
+	if n != len(right) {
 		return newErrorFormulaArg(formulaErrorNA, formulaErrorNA)
 	}
-	number, val, cnt := 0.0, 0.0, 0.0
-	for token := argsList.Front(); token != nil; token = token.Next() {
-		arg := token.Value.(formulaArg)
+	l1, l2 := list.New(), list.New()
+	l1.PushBack(array1)
+	l2.PushBack(array2)
+	result, skip := 0.0, 0
+	mean1, mean2 := fn.AVERAGE(l1), fn.AVERAGE(l2)
+	for i := 0; i < n; i++ {
+		arg1 := left[i].ToNumber()
+		arg2 := right[i].ToNumber()
+		if arg1.Type == ArgError || arg2.Type == ArgError {
+			skip++
+			continue
+		}
+		result += (arg1.Number - mean1.Number) * (arg2.Number - mean2.Number)
+	}
+	return newNumberFormulaArg(result / float64(n-skip))
+}
+
+// COVARIANCEdotP function calculates the population covariance of two supplied
+// sets of values. The syntax of the function is:
+//
+//    COVARIANCE.P(array1,array2)
+//
+func (fn *formulaFuncs) COVARIANCEdotP(argsList *list.List) formulaArg {
+	if argsList.Len() != 2 {
+		return newErrorFormulaArg(formulaErrorVALUE, "COVARIANCE.P requires 2 arguments")
+	}
+	return fn.COVAR(argsList)
+}
+
+// calcStringCountSum is part of the implementation countSum.
+func calcStringCountSum(countText bool, count, sum float64, num, arg formulaArg) (float64, float64) {
+	if countText && num.Type == ArgError && arg.String != "" {
+		count++
+	}
+	if num.Type == ArgNumber {
+		sum += num.Number
+		count++
+	}
+	return count, sum
+}
+
+// countSum get count and sum for a formula arguments array.
+func (fn *formulaFuncs) countSum(countText bool, args []formulaArg) (count, sum float64) {
+	for _, arg := range args {
 		switch arg.Type {
+		case ArgNumber:
+			if countText || !arg.Boolean {
+				sum += arg.Number
+				count++
+			}
 		case ArgString:
-			num := arg.ToNumber()
-			if num.Type != ArgNumber {
+			if !countText && (arg.Value() == "TRUE" || arg.Value() == "FALSE") {
 				continue
+			} else if countText && (arg.Value() == "TRUE" || arg.Value() == "FALSE") {
+				num := arg.ToBool()
+				if num.Type == ArgNumber {
+					count++
+					sum += num.Number
+					continue
+				}
 			}
-			number = num.Number
-		case ArgNumber:
-			number = arg.Number
-		}
-		if number <= 0 {
-			return newErrorFormulaArg(formulaErrorNA, formulaErrorNA)
+			num := arg.ToNumber()
+			count, sum = calcStringCountSum(countText, count, sum, num, arg)
+		case ArgList, ArgMatrix:
+			cnt, summary := fn.countSum(countText, arg.ToList())
+			sum += summary
+			count += cnt
 		}
-		val += (1 / number)
-		cnt++
 	}
-	return newNumberFormulaArg(1 / (val / cnt))
+	return
 }
 
-// KURT function calculates the kurtosis of a supplied set of values. The
-// syntax of the function is:
+// CORREL function calculates the Pearson Product-Moment Correlation
+// Coefficient for two sets of values. The syntax of the function is:
 //
-//    KURT(number1,[number2],...)
+//    CORREL(array1,array2)
 //
-func (fn *formulaFuncs) KURT(argsList *list.List) formulaArg {
-	if argsList.Len() < 1 {
-		return newErrorFormulaArg(formulaErrorVALUE, "KURT requires at least 1 argument")
+func (fn *formulaFuncs) CORREL(argsList *list.List) formulaArg {
+	if argsList.Len() != 2 {
+		return newErrorFormulaArg(formulaErrorVALUE, "CORREL requires 2 arguments")
 	}
-	mean, stdev := fn.AVERAGE(argsList), fn.STDEV(argsList)
-	if stdev.Number > 0 {
-		count, summer := 0.0, 0.0
-		for arg := argsList.Front(); arg != nil; arg = arg.Next() {
-			token := arg.Value.(formulaArg)
-			switch token.Type {
-			case ArgString, ArgNumber:
-				num := token.ToNumber()
-				if num.Type == ArgError {
-					continue
-				}
-				summer += math.Pow((num.Number-mean.Number)/stdev.Number, 4)
+	array1 := argsList.Front().Value.(formulaArg)
+	array2 := argsList.Back().Value.(formulaArg)
+	left, right := array1.ToList(), array2.ToList()
+	n := len(left)
+	if n != len(right) {
+		return newErrorFormulaArg(formulaErrorNA, formulaErrorNA)
+	}
+	l1, l2, l3 := list.New(), list.New(), list.New()
+	for i := 0; i < n; i++ {
+		if lhs, rhs := left[i].ToNumber(), right[i].ToNumber(); lhs.Number != 0 && rhs.Number != 0 {
+			l1.PushBack(lhs)
+			l2.PushBack(rhs)
+		}
+	}
+	stdev1, stdev2 := fn.STDEV(l1), fn.STDEV(l2)
+	if stdev1.Number == 0 || stdev2.Number == 0 {
+		return newErrorFormulaArg(formulaErrorDIV, formulaErrorDIV)
+	}
+	mean1, mean2, skip := fn.AVERAGE(l1), fn.AVERAGE(l2), 0
+	for i := 0; i < n; i++ {
+		lhs, rhs := left[i].ToNumber(), right[i].ToNumber()
+		if lhs.Number == 0 || rhs.Number == 0 {
+			skip++
+			continue
+		}
+		l3.PushBack(newNumberFormulaArg((lhs.Number - mean1.Number) * (rhs.Number - mean2.Number)))
+	}
+	return newNumberFormulaArg(fn.SUM(l3).Number / float64(n-skip-1) / stdev1.Number / stdev2.Number)
+}
+
+// COUNT function returns the count of numeric values in a supplied set of
+// cells or values. This count includes both numbers and dates. The syntax of
+// the function is:
+//
+//    COUNT(value1,[value2],...)
+//
+func (fn *formulaFuncs) COUNT(argsList *list.List) formulaArg {
+	var count int
+	for token := argsList.Front(); token != nil; token = token.Next() {
+		arg := token.Value.(formulaArg)
+		switch arg.Type {
+		case ArgString:
+			if arg.ToNumber().Type != ArgError {
 				count++
-			case ArgList, ArgMatrix:
-				for _, row := range token.ToList() {
-					if row.Type == ArgNumber || row.Type == ArgString {
-						num := row.ToNumber()
-						if num.Type == ArgError {
-							continue
-						}
-						summer += math.Pow((num.Number-mean.Number)/stdev.Number, 4)
+			}
+		case ArgNumber:
+			count++
+		case ArgMatrix:
+			for _, row := range arg.Matrix {
+				for _, value := range row {
+					if value.ToNumber().Type != ArgError {
 						count++
 					}
 				}
 			}
 		}
-		if count > 3 {
-			return newNumberFormulaArg(summer*(count*(count+1)/((count-1)*(count-2)*(count-3))) - (3 * math.Pow(count-1, 2) / ((count - 2) * (count - 3))))
-		}
 	}
-	return newErrorFormulaArg(formulaErrorDIV, formulaErrorDIV)
+	return newNumberFormulaArg(float64(count))
 }
 
-// EXPONdotDIST function returns the value of the exponential distribution for
-// a give value of x. The user can specify whether the probability density
-// function or the cumulative distribution function is used. The syntax of the
-// Expondist function is:
+// COUNTA function returns the number of non-blanks within a supplied set of
+// cells or values. The syntax of the function is:
 //
-//    EXPON.DIST(x,lambda,cumulative)
+//    COUNTA(value1,[value2],...)
 //
-func (fn *formulaFuncs) EXPONdotDIST(argsList *list.List) formulaArg {
-	if argsList.Len() != 3 {
-		return newErrorFormulaArg(formulaErrorVALUE, "EXPON.DIST requires 3 arguments")
+func (fn *formulaFuncs) COUNTA(argsList *list.List) formulaArg {
+	var count int
+	for token := argsList.Front(); token != nil; token = token.Next() {
+		arg := token.Value.(formulaArg)
+		switch arg.Type {
+		case ArgString:
+			if arg.String != "" {
+				count++
+			}
+		case ArgNumber:
+			count++
+		case ArgMatrix:
+			for _, row := range arg.ToList() {
+				switch row.Type {
+				case ArgString:
+					if row.String != "" {
+						count++
+					}
+				case ArgNumber:
+					count++
+				}
+			}
+		}
 	}
-	return fn.EXPONDIST(argsList)
+	return newNumberFormulaArg(float64(count))
 }
 
-// EXPONDIST function returns the value of the exponential distribution for a
-// give value of x. The user can specify whether the probability density
-// function or the cumulative distribution function is used. The syntax of the
-// Expondist function is:
+// COUNTBLANK function returns the number of blank cells in a supplied range.
+// The syntax of the function is:
 //
-//    EXPONDIST(x,lambda,cumulative)
+//    COUNTBLANK(range)
 //
-func (fn *formulaFuncs) EXPONDIST(argsList *list.List) formulaArg {
-	if argsList.Len() != 3 {
-		return newErrorFormulaArg(formulaErrorVALUE, "EXPONDIST requires 3 arguments")
-	}
-	var x, lambda, cumulative formulaArg
-	if x = argsList.Front().Value.(formulaArg).ToNumber(); x.Type != ArgNumber {
-		return x
-	}
-	if lambda = argsList.Front().Next().Value.(formulaArg).ToNumber(); lambda.Type != ArgNumber {
-		return lambda
+func (fn *formulaFuncs) COUNTBLANK(argsList *list.List) formulaArg {
+	if argsList.Len() != 1 {
+		return newErrorFormulaArg(formulaErrorVALUE, "COUNTBLANK requires 1 argument")
 	}
-	if cumulative = argsList.Back().Value.(formulaArg).ToBool(); cumulative.Type == ArgError {
-		return cumulative
-	}
-	if x.Number < 0 {
-		return newErrorFormulaArg(formulaErrorNUM, formulaErrorNUM)
-	}
-	if lambda.Number <= 0 {
-		return newErrorFormulaArg(formulaErrorNUM, formulaErrorNUM)
-	}
-	if cumulative.Number == 1 {
-		return newNumberFormulaArg(1 - math.Exp(-lambda.Number*x.Number))
+	var count float64
+	for _, cell := range argsList.Front().Value.(formulaArg).ToList() {
+		if cell.Value() == "" {
+			count++
+		}
 	}
-	return newNumberFormulaArg(lambda.Number * math.Exp(-lambda.Number*x.Number))
+	return newNumberFormulaArg(count)
 }
 
-// d1mach returns double precision real machine constants.
-func d1mach(i int) float64 {
-	arr := []float64{
-		2.2250738585072014e-308,
-		1.7976931348623158e+308,
-		1.1102230246251565e-16,
-		2.2204460492503131e-16,
-		0.301029995663981195,
+// COUNTIF function returns the number of cells within a supplied range, that
+// satisfy a given criteria. The syntax of the function is:
+//
+//    COUNTIF(range,criteria)
+//
+func (fn *formulaFuncs) COUNTIF(argsList *list.List) formulaArg {
+	if argsList.Len() != 2 {
+		return newErrorFormulaArg(formulaErrorVALUE, "COUNTIF requires 2 arguments")
 	}
-	if i > len(arr) {
-		return 0
+	var (
+		criteria = formulaCriteriaParser(argsList.Front().Next().Value.(formulaArg).String)
+		count    float64
+	)
+	for _, cell := range argsList.Front().Value.(formulaArg).ToList() {
+		if ok, _ := formulaCriteriaEval(cell.Value(), criteria); ok {
+			count++
+		}
 	}
-	return arr[i-1]
+	return newNumberFormulaArg(count)
 }
 
-// chebyshevInit determines the number of terms for the double precision
-// orthogonal series "dos" needed to insure the error is no larger
-// than "eta". Ordinarily eta will be chosen to be one-tenth machine
-// precision.
-func chebyshevInit(nos int, eta float64, dos []float64) int {
-	i, e := 0, 0.0
-	if nos < 1 {
-		return 0
-	}
-	for ii := 1; ii <= nos; ii++ {
-		i = nos - ii
-		e += math.Abs(dos[i])
-		if e > eta {
-			return i
+// formulaIfsMatch function returns cells reference array which match criteria.
+func formulaIfsMatch(args []formulaArg) (cellRefs []cellRef) {
+	for i := 0; i < len(args)-1; i += 2 {
+		match := []cellRef{}
+		matrix, criteria := args[i].Matrix, formulaCriteriaParser(args[i+1].Value())
+		if i == 0 {
+			for rowIdx, row := range matrix {
+				for colIdx, col := range row {
+					if ok, _ := formulaCriteriaEval(col.Value(), criteria); ok {
+						match = append(match, cellRef{Col: colIdx, Row: rowIdx})
+					}
+				}
+			}
+		} else {
+			for _, ref := range cellRefs {
+				value := matrix[ref.Row][ref.Col]
+				if ok, _ := formulaCriteriaEval(value.Value(), criteria); ok {
+					match = append(match, ref)
+				}
+			}
+		}
+		if len(match) == 0 {
+			return
 		}
+		cellRefs = match[:]
 	}
-	return i
+	return
 }
 
-// chebyshevEval evaluates the n-term Chebyshev series "a" at "x".
-func chebyshevEval(n int, x float64, a []float64) float64 {
-	if n < 1 || n > 1000 || x < -1.1 || x > 1.1 {
-		return math.NaN()
-	}
-	twox, b0, b1, b2 := x*2, 0.0, 0.0, 0.0
-	for i := 1; i <= n; i++ {
-		b2 = b1
-		b1 = b0
-		b0 = twox*b1 - b2 + a[n-i]
+// COUNTIFS function returns the number of rows within a table, that satisfy a
+// set of given criteria. The syntax of the function is:
+//
+//    COUNTIFS(criteria_range1,criteria1,[criteria_range2,criteria2],...)
+//
+func (fn *formulaFuncs) COUNTIFS(argsList *list.List) formulaArg {
+	if argsList.Len() < 2 {
+		return newErrorFormulaArg(formulaErrorVALUE, "COUNTIFS requires at least 2 arguments")
 	}
-	return (b0 - b2) * 0.5
-}
-
-// lgammacor is an implementation for the log(gamma) correction.
-func lgammacor(x float64) float64 {
-	algmcs := []float64{
-		0.1666389480451863247205729650822, -0.1384948176067563840732986059135e-4,
-		0.9810825646924729426157171547487e-8, -0.1809129475572494194263306266719e-10,
-		0.6221098041892605227126015543416e-13, -0.3399615005417721944303330599666e-15,
-		0.2683181998482698748957538846666e-17, -0.2868042435334643284144622399999e-19,
-		0.3962837061046434803679306666666e-21, -0.6831888753985766870111999999999e-23,
-		0.1429227355942498147573333333333e-24, -0.3547598158101070547199999999999e-26,
-		0.1025680058010470912000000000000e-27, -0.3401102254316748799999999999999e-29,
-		0.1276642195630062933333333333333e-30,
+	if argsList.Len()%2 != 0 {
+		return newErrorFormulaArg(formulaErrorNA, formulaErrorNA)
 	}
-	nalgm := chebyshevInit(15, d1mach(3), algmcs)
-	xbig := 1.0 / math.Sqrt(d1mach(3))
-	xmax := math.Exp(math.Min(math.Log(d1mach(2)/12.0), -math.Log(12.0*d1mach(1))))
-	if x < 10.0 {
-		return math.NaN()
-	} else if x >= xmax {
-		return 4.930380657631324e-32
-	} else if x < xbig {
-		tmp := 10.0 / x
-		return chebyshevEval(nalgm, tmp*tmp*2.0-1.0, algmcs) / x
+	args := []formulaArg{}
+	for arg := argsList.Front(); arg != nil; arg = arg.Next() {
+		args = append(args, arg.Value.(formulaArg))
 	}
-	return 1.0 / (x * 12.0)
+	return newNumberFormulaArg(float64(len(formulaIfsMatch(args))))
 }
 
-// logrelerr compute the relative error logarithm.
-func logrelerr(x float64) float64 {
-	alnrcs := []float64{
-		0.10378693562743769800686267719098e+1, -0.13364301504908918098766041553133,
-		0.19408249135520563357926199374750e-1, -0.30107551127535777690376537776592e-2,
-		0.48694614797154850090456366509137e-3, -0.81054881893175356066809943008622e-4,
-		0.13778847799559524782938251496059e-4, -0.23802210894358970251369992914935e-5,
-		0.41640416213865183476391859901989e-6, -0.73595828378075994984266837031998e-7,
-		0.13117611876241674949152294345011e-7, -0.23546709317742425136696092330175e-8,
-		0.42522773276034997775638052962567e-9, -0.77190894134840796826108107493300e-10,
-		0.14075746481359069909215356472191e-10, -0.25769072058024680627537078627584e-11,
-		0.47342406666294421849154395005938e-12, -0.87249012674742641745301263292675e-13,
-		0.16124614902740551465739833119115e-13, -0.29875652015665773006710792416815e-14,
-		0.55480701209082887983041321697279e-15, -0.10324619158271569595141333961932e-15,
-		0.19250239203049851177878503244868e-16, -0.35955073465265150011189707844266e-17,
-		0.67264542537876857892194574226773e-18, -0.12602624168735219252082425637546e-18,
-		0.23644884408606210044916158955519e-19, -0.44419377050807936898878389179733e-20,
-		0.83546594464034259016241293994666e-21, -0.15731559416479562574899253521066e-21,
-		0.29653128740247422686154369706666e-22, -0.55949583481815947292156013226666e-23,
-		0.10566354268835681048187284138666e-23, -0.19972483680670204548314999466666e-24,
-		0.37782977818839361421049855999999e-25, -0.71531586889081740345038165333333e-26,
-		0.13552488463674213646502024533333e-26, -0.25694673048487567430079829333333e-27,
-		0.48747756066216949076459519999999e-28, -0.92542112530849715321132373333333e-29,
-		0.17578597841760239233269760000000e-29, -0.33410026677731010351377066666666e-30,
-		0.63533936180236187354180266666666e-31,
+// DEVSQ function calculates the sum of the squared deviations from the sample
+// mean. The syntax of the function is:
+//
+//    DEVSQ(number1,[number2],...)
+//
+func (fn *formulaFuncs) DEVSQ(argsList *list.List) formulaArg {
+	if argsList.Len() < 1 {
+		return newErrorFormulaArg(formulaErrorVALUE, "DEVSQ requires at least 1 numeric argument")
 	}
-	nlnrel := chebyshevInit(43, 0.1*d1mach(3), alnrcs)
-	if x <= -1 {
-		return math.NaN()
+	avg, count, result := fn.AVERAGE(argsList), -1, 0.0
+	for arg := argsList.Front(); arg != nil; arg = arg.Next() {
+		for _, number := range arg.Value.(formulaArg).ToList() {
+			num := number.ToNumber()
+			if num.Type != ArgNumber {
+				continue
+			}
+			count++
+			if count == 0 {
+				result = math.Pow(num.Number-avg.Number, 2)
+				continue
+			}
+			result += math.Pow(num.Number-avg.Number, 2)
+		}
 	}
-	if math.Abs(x) <= 0.375 {
-		return x * (1.0 - x*chebyshevEval(nlnrel, x/0.375, alnrcs))
+	if count == -1 {
+		return newErrorFormulaArg(formulaErrorNA, formulaErrorNA)
 	}
-	return math.Log(x + 1.0)
+	return newNumberFormulaArg(result)
 }
 
-// logBeta is an implementation for the log of the beta distribution
-// function.
-func logBeta(a, b float64) float64 {
-	corr, p, q := 0.0, a, a
-	if b < p {
-		p = b
+// FISHER function calculates the Fisher Transformation for a supplied value.
+// The syntax of the function is:
+//
+//    FISHER(x)
+//
+func (fn *formulaFuncs) FISHER(argsList *list.List) formulaArg {
+	if argsList.Len() != 1 {
+		return newErrorFormulaArg(formulaErrorVALUE, "FISHER requires 1 numeric argument")
 	}
-	if b > q {
-		q = b
+	token := argsList.Front().Value.(formulaArg)
+	switch token.Type {
+	case ArgString:
+		arg := token.ToNumber()
+		if arg.Type == ArgNumber {
+			if arg.Number <= -1 || arg.Number >= 1 {
+				return newErrorFormulaArg(formulaErrorNA, formulaErrorNA)
+			}
+			return newNumberFormulaArg(0.5 * math.Log((1+arg.Number)/(1-arg.Number)))
+		}
+	case ArgNumber:
+		if token.Number <= -1 || token.Number >= 1 {
+			return newErrorFormulaArg(formulaErrorNA, formulaErrorNA)
+		}
+		return newNumberFormulaArg(0.5 * math.Log((1+token.Number)/(1-token.Number)))
 	}
-	if p < 0 {
-		return math.NaN()
+	return newErrorFormulaArg(formulaErrorVALUE, "FISHER requires 1 numeric argument")
+}
+
+// FISHERINV function calculates the inverse of the Fisher Transformation and
+// returns a value between -1 and +1. The syntax of the function is:
+//
+//    FISHERINV(y)
+//
+func (fn *formulaFuncs) FISHERINV(argsList *list.List) formulaArg {
+	if argsList.Len() != 1 {
+		return newErrorFormulaArg(formulaErrorVALUE, "FISHERINV requires 1 numeric argument")
 	}
-	if p == 0 {
-		return math.MaxFloat64
+	token := argsList.Front().Value.(formulaArg)
+	switch token.Type {
+	case ArgString:
+		arg := token.ToNumber()
+		if arg.Type == ArgNumber {
+			return newNumberFormulaArg((math.Exp(2*arg.Number) - 1) / (math.Exp(2*arg.Number) + 1))
+		}
+	case ArgNumber:
+		return newNumberFormulaArg((math.Exp(2*token.Number) - 1) / (math.Exp(2*token.Number) + 1))
 	}
-	if p >= 10.0 {
-		corr = lgammacor(p) + lgammacor(q) - lgammacor(p+q)
-		return math.Log(q)*-0.5 + 0.918938533204672741780329736406 + corr + (p-0.5)*math.Log(p/(p+q)) + q*logrelerr(-p/(p+q))
+	return newErrorFormulaArg(formulaErrorVALUE, "FISHERINV requires 1 numeric argument")
+}
+
+// GAMMA function returns the value of the Gamma Function, Γ(n), for a
+// specified number, n. The syntax of the function is:
+//
+//    GAMMA(number)
+//
+func (fn *formulaFuncs) GAMMA(argsList *list.List) formulaArg {
+	if argsList.Len() != 1 {
+		return newErrorFormulaArg(formulaErrorVALUE, "GAMMA requires 1 numeric argument")
 	}
-	if q >= 10 {
-		corr = lgammacor(q) - lgammacor(p+q)
-		val, _ := math.Lgamma(p)
-		return val + corr + p - p*math.Log(p+q) + (q-0.5)*logrelerr(-p/(p+q))
+	token := argsList.Front().Value.(formulaArg)
+	switch token.Type {
+	case ArgString:
+		arg := token.ToNumber()
+		if arg.Type == ArgNumber {
+			if arg.Number <= 0 {
+				return newErrorFormulaArg(formulaErrorNA, formulaErrorNA)
+			}
+			return newNumberFormulaArg(math.Gamma(arg.Number))
+		}
+	case ArgNumber:
+		if token.Number <= 0 {
+			return newErrorFormulaArg(formulaErrorNA, formulaErrorNA)
+		}
+		return newNumberFormulaArg(math.Gamma(token.Number))
 	}
-	return math.Log(math.Gamma(p) * (math.Gamma(q) / math.Gamma(p+q)))
+	return newErrorFormulaArg(formulaErrorVALUE, "GAMMA requires 1 numeric argument")
 }
 
-// pbetaRaw is a part of pbeta for the beta distribution.
-func pbetaRaw(alnsml, ans, eps, p, pin, q, sml, x, y float64) float64 {
-	if q > 1.0 {
-		xb := p*math.Log(y) + q*math.Log(1.0-y) - logBeta(p, q) - math.Log(q)
-		ib := int(math.Max(xb/alnsml, 0.0))
-		term := math.Exp(xb - float64(ib)*alnsml)
-		c := 1.0 / (1.0 - y)
-		p1 := q * c / (p + q - 1.0)
-		finsum := 0.0
-		n := int(q)
-		if q == float64(n) {
-			n = n - 1
-		}
-		for i := 1; i <= n; i++ {
-			if p1 <= 1 && term/eps <= finsum {
-				break
-			}
-			xi := float64(i)
-			term = (q - xi + 1.0) * c * term / (p + q - xi)
-			if term > 1.0 {
-				ib = ib - 1
-				term = term * sml
-			}
-			if ib == 0 {
-				finsum = finsum + term
+// GAMMALN function returns the natural logarithm of the Gamma Function, Γ
+// (n). The syntax of the function is:
+//
+//    GAMMALN(x)
+//
+func (fn *formulaFuncs) GAMMALN(argsList *list.List) formulaArg {
+	if argsList.Len() != 1 {
+		return newErrorFormulaArg(formulaErrorVALUE, "GAMMALN requires 1 numeric argument")
+	}
+	token := argsList.Front().Value.(formulaArg)
+	switch token.Type {
+	case ArgString:
+		arg := token.ToNumber()
+		if arg.Type == ArgNumber {
+			if arg.Number <= 0 {
+				return newErrorFormulaArg(formulaErrorNA, formulaErrorNA)
 			}
+			return newNumberFormulaArg(math.Log(math.Gamma(arg.Number)))
+		}
+	case ArgNumber:
+		if token.Number <= 0 {
+			return newErrorFormulaArg(formulaErrorNA, formulaErrorNA)
 		}
-		ans = ans + finsum
+		return newNumberFormulaArg(math.Log(math.Gamma(token.Number)))
 	}
-	if y != x || p != pin {
-		ans = 1.0 - ans
+	return newErrorFormulaArg(formulaErrorVALUE, "GAMMALN requires 1 numeric argument")
+}
+
+// GEOMEAN function calculates the geometric mean of a supplied set of values.
+// The syntax of the function is:
+//
+//    GEOMEAN(number1,[number2],...)
+//
+func (fn *formulaFuncs) GEOMEAN(argsList *list.List) formulaArg {
+	if argsList.Len() < 1 {
+		return newErrorFormulaArg(formulaErrorVALUE, "GEOMEAN requires at least 1 numeric argument")
 	}
-	ans = math.Max(math.Min(ans, 1.0), 0.0)
-	return ans
+	product := fn.PRODUCT(argsList)
+	if product.Type != ArgNumber {
+		return product
+	}
+	count := fn.COUNT(argsList)
+	min := fn.MIN(argsList)
+	if product.Number > 0 && min.Number > 0 {
+		return newNumberFormulaArg(math.Pow(product.Number, (1 / count.Number)))
+	}
+	return newErrorFormulaArg(formulaErrorNUM, formulaErrorNUM)
 }
 
-// pbeta returns distribution function of the beta distribution.
-func pbeta(x, pin, qin float64) (ans float64) {
-	eps := d1mach(3)
-	alneps := math.Log(eps)
-	sml := d1mach(1)
-	alnsml := math.Log(sml)
-	y := x
-	p := pin
-	q := qin
-	if p/(p+q) < x {
-		y = 1.0 - y
-		p = qin
-		q = pin
+// HARMEAN function calculates the harmonic mean of a supplied set of values.
+// The syntax of the function is:
+//
+//    HARMEAN(number1,[number2],...)
+//
+func (fn *formulaFuncs) HARMEAN(argsList *list.List) formulaArg {
+	if argsList.Len() < 1 {
+		return newErrorFormulaArg(formulaErrorVALUE, "HARMEAN requires at least 1 argument")
 	}
-	if (p+q)*y/(p+1.0) < eps {
-		xb := p*math.Log(math.Max(y, sml)) - math.Log(p) - logBeta(p, q)
-		if xb > alnsml && y != 0.0 {
-			ans = math.Exp(xb)
-		}
-		if y != x || p != pin {
-			ans = 1.0 - ans
-		}
-	} else {
-		ps := q - math.Floor(q)
-		if ps == 0.0 {
-			ps = 1.0
-		}
-		xb := p*math.Log(y) - logBeta(ps, p) - math.Log(p)
-		if xb >= alnsml {
-			ans = math.Exp(xb)
-			term := ans * p
-			if ps != 1.0 {
-				n := int(math.Max(alneps/math.Log(y), 4.0))
-				for i := 1; i <= n; i++ {
-					xi := float64(i)
-					term = term * (xi - ps) * y / xi
-					ans = ans + term/(p+xi)
-				}
+	if min := fn.MIN(argsList); min.Number < 0 {
+		return newErrorFormulaArg(formulaErrorNA, formulaErrorNA)
+	}
+	number, val, cnt := 0.0, 0.0, 0.0
+	for token := argsList.Front(); token != nil; token = token.Next() {
+		arg := token.Value.(formulaArg)
+		switch arg.Type {
+		case ArgString:
+			num := arg.ToNumber()
+			if num.Type != ArgNumber {
+				continue
 			}
+			number = num.Number
+		case ArgNumber:
+			number = arg.Number
 		}
-		ans = pbetaRaw(alnsml, ans, eps, p, pin, q, sml, x, y)
+		if number <= 0 {
+			return newErrorFormulaArg(formulaErrorNA, formulaErrorNA)
+		}
+		val += (1 / number)
+		cnt++
 	}
-	return ans
+	return newNumberFormulaArg(1 / (val / cnt))
 }
 
-// betainvProbIterator is a part of betainv for the inverse of the beta function.
-func betainvProbIterator(alpha1, alpha3, beta1, beta2, beta3, logbeta, lower, maxCumulative, prob1, prob2, upper float64, needSwap bool) float64 {
-	var i, j, prev, prop4 float64
-	j = 1
-	for prob := 0; prob < 1000; prob++ {
-		prop3 := pbeta(beta3, alpha1, beta1)
-		prop3 = (prop3 - prob1) * math.Exp(logbeta+prob2*math.Log(beta3)+beta2*math.Log(1.0-beta3))
-		if prop3*prop4 <= 0 {
-			prev = math.Max(math.Abs(j), maxCumulative)
-		}
-		h := 1.0
-		for iteratorCount := 0; iteratorCount < 1000; iteratorCount++ {
-			j = h * prop3
-			if math.Abs(j) < prev {
-				i = beta3 - j
-				if i >= 0 && i <= 1.0 {
-					if prev <= alpha3 {
-						return beta3
-					}
-					if math.Abs(prop3) <= alpha3 {
-						return beta3
-					}
-					if i != 0 && i != 1.0 {
-						break
+// KURT function calculates the kurtosis of a supplied set of values. The
+// syntax of the function is:
+//
+//    KURT(number1,[number2],...)
+//
+func (fn *formulaFuncs) KURT(argsList *list.List) formulaArg {
+	if argsList.Len() < 1 {
+		return newErrorFormulaArg(formulaErrorVALUE, "KURT requires at least 1 argument")
+	}
+	mean, stdev := fn.AVERAGE(argsList), fn.STDEV(argsList)
+	if stdev.Number > 0 {
+		count, summer := 0.0, 0.0
+		for arg := argsList.Front(); arg != nil; arg = arg.Next() {
+			token := arg.Value.(formulaArg)
+			switch token.Type {
+			case ArgString, ArgNumber:
+				num := token.ToNumber()
+				if num.Type == ArgError {
+					continue
+				}
+				summer += math.Pow((num.Number-mean.Number)/stdev.Number, 4)
+				count++
+			case ArgList, ArgMatrix:
+				for _, row := range token.ToList() {
+					if row.Type == ArgNumber || row.Type == ArgString {
+						num := row.ToNumber()
+						if num.Type == ArgError {
+							continue
+						}
+						summer += math.Pow((num.Number-mean.Number)/stdev.Number, 4)
+						count++
 					}
 				}
 			}
-			h /= 3.0
 		}
-		if i == beta3 {
-			return beta3
+		if count > 3 {
+			return newNumberFormulaArg(summer*(count*(count+1)/((count-1)*(count-2)*(count-3))) - (3 * math.Pow(count-1, 2) / ((count - 2) * (count - 3))))
 		}
-		beta3, prop4 = i, prop3
 	}
-	return beta3
+	return newErrorFormulaArg(formulaErrorDIV, formulaErrorDIV)
 }
 
-// betainv is an implementation for the quantile of the beta distribution.
-func betainv(probability, alpha, beta, lower, upper float64) float64 {
-	minCumulative, maxCumulative := 1.0e-300, 3.0e-308
-	lowerBound, upperBound := maxCumulative, 1.0-2.22e-16
-	needSwap := false
-	var alpha1, alpha2, beta1, beta2, beta3, prob1, x, y float64
-	if probability <= 0.5 {
-		prob1, alpha1, beta1 = probability, alpha, beta
-	} else {
-		prob1, alpha1, beta1, needSwap = 1.0-probability, beta, alpha, true
+// EXPONdotDIST function returns the value of the exponential distribution for
+// a give value of x. The user can specify whether the probability density
+// function or the cumulative distribution function is used. The syntax of the
+// Expondist function is:
+//
+//    EXPON.DIST(x,lambda,cumulative)
+//
+func (fn *formulaFuncs) EXPONdotDIST(argsList *list.List) formulaArg {
+	if argsList.Len() != 3 {
+		return newErrorFormulaArg(formulaErrorVALUE, "EXPON.DIST requires 3 arguments")
 	}
-	logbeta := logBeta(alpha, beta)
-	prob2 := math.Sqrt(-math.Log(prob1 * prob1))
-	prob3 := prob2 - (prob2*0.27061+2.3075)/(prob2*(prob2*0.04481+0.99229)+1)
-	if alpha1 > 1 && beta1 > 1 {
-		alpha2, beta2, prob2 = 1/(alpha1+alpha1-1), 1/(beta1+beta1-1), (prob3*prob3-3)/6
-		x = 2 / (alpha2 + beta2)
-		y = prob3*math.Sqrt(x+prob2)/x - (beta2-alpha2)*(prob2+5/6.0-2/(x*3))
-		beta3 = alpha1 / (alpha1 + beta1*math.Exp(y+y))
-	} else {
-		beta2, prob2 = 1/(beta1*9), beta1+beta1
-		beta2 = prob2 * math.Pow(1-beta2+prob3*math.Sqrt(beta2), 3)
-		if beta2 <= 0 {
-			beta3 = 1 - math.Exp((math.Log((1-prob1)*beta1)+logbeta)/beta1)
-		} else {
-			beta2 = (prob2 + alpha1*4 - 2) / beta2
-			if beta2 <= 1 {
-				beta3 = math.Exp((logbeta + math.Log(alpha1*prob1)) / alpha1)
-			} else {
-				beta3 = 1 - 2/(beta2+1)
-			}
-		}
+	return fn.EXPONDIST(argsList)
+}
+
+// EXPONDIST function returns the value of the exponential distribution for a
+// give value of x. The user can specify whether the probability density
+// function or the cumulative distribution function is used. The syntax of the
+// Expondist function is:
+//
+//    EXPONDIST(x,lambda,cumulative)
+//
+func (fn *formulaFuncs) EXPONDIST(argsList *list.List) formulaArg {
+	if argsList.Len() != 3 {
+		return newErrorFormulaArg(formulaErrorVALUE, "EXPONDIST requires 3 arguments")
 	}
-	beta2, prob2 = 1-beta1, 1-alpha1
-	if beta3 < lowerBound {
-		beta3 = lowerBound
-	} else if beta3 > upperBound {
-		beta3 = upperBound
+	var x, lambda, cumulative formulaArg
+	if x = argsList.Front().Value.(formulaArg).ToNumber(); x.Type != ArgNumber {
+		return x
 	}
-	alpha3 := math.Max(minCumulative, math.Pow(10.0, -13.0-2.5/(alpha1*alpha1)-0.5/(prob1*prob1)))
-	beta3 = betainvProbIterator(alpha1, alpha3, beta1, beta2, beta3, logbeta, lower, maxCumulative, prob1, prob2, upper, needSwap)
-	if needSwap {
-		beta3 = 1.0 - beta3
+	if lambda = argsList.Front().Next().Value.(formulaArg).ToNumber(); lambda.Type != ArgNumber {
+		return lambda
 	}
-	return (upper-lower)*beta3 + lower
+	if cumulative = argsList.Back().Value.(formulaArg).ToBool(); cumulative.Type == ArgError {
+		return cumulative
+	}
+	if x.Number < 0 {
+		return newErrorFormulaArg(formulaErrorNUM, formulaErrorNUM)
+	}
+	if lambda.Number <= 0 {
+		return newErrorFormulaArg(formulaErrorNUM, formulaErrorNUM)
+	}
+	if cumulative.Number == 1 {
+		return newNumberFormulaArg(1 - math.Exp(-lambda.Number*x.Number))
+	}
+	return newNumberFormulaArg(lambda.Number * math.Exp(-lambda.Number*x.Number))
 }
 
-// FINV function calculates the inverse of the (right-tailed) F Probability
-// Distribution for a supplied probability. The syntax of the function is:
-//
-//    FINV(probability,deg_freedom1,deg_freedom2)
-//
-func (fn *formulaFuncs) FINV(argsList *list.List) formulaArg {
+// finv is an implementation of the formula functions F.INV.RT and FINV.
+func (fn *formulaFuncs) finv(name string, argsList *list.List) formulaArg {
 	if argsList.Len() != 3 {
-		return newErrorFormulaArg(formulaErrorVALUE, "FINV requires 3 arguments")
+		return newErrorFormulaArg(formulaErrorVALUE, fmt.Sprintf("%s requires 3 arguments", name))
 	}
 	var probability, d1, d2 formulaArg
 	if probability = argsList.Front().Value.(formulaArg).ToNumber(); probability.Type != ArgNumber {
@@ -6167,7 +6232,26 @@ func (fn *formulaFuncs) FINV(argsList *list.List) formulaArg {
 	if d2.Number < 1 || d2.Number >= math.Pow10(10) {
 		return newErrorFormulaArg(formulaErrorNUM, formulaErrorNUM)
 	}
-	return newNumberFormulaArg((1/betainv(1.0-(1.0-probability.Number), d2.Number/2, d1.Number/2, 0, 1) - 1.0) * (d2.Number / d1.Number))
+	return newNumberFormulaArg((1/calcBetainv(1.0-(1.0-probability.Number), d2.Number/2, d1.Number/2, 0, 1) - 1.0) * (d2.Number / d1.Number))
+}
+
+// FdotINVdotRT function calculates the inverse of the (right-tailed) F
+// Probability Distribution for a supplied probability. The syntax of the
+// function is:
+//
+//    F.INV.RT(probability,deg_freedom1,deg_freedom2)
+//
+func (fn *formulaFuncs) FdotINVdotRT(argsList *list.List) formulaArg {
+	return fn.finv("F.INV.RT", argsList)
+}
+
+// FINV function calculates the inverse of the (right-tailed) F Probability
+// Distribution for a supplied probability. The syntax of the function is:
+//
+//    FINV(probability,deg_freedom1,deg_freedom2)
+//
+func (fn *formulaFuncs) FINV(argsList *list.List) formulaArg {
+	return fn.finv("FINV", argsList)
 }
 
 // NORMdotDIST function calculates the Normal Probability Density Function or
diff --git a/calc_test.go b/calc_test.go
index 80ba8ef774fabcc72532c95bd3afb2e48c300958..63e5984efb21c66f18e0015d5ec52060ce0d0ddc 100644
--- a/calc_test.go
+++ b/calc_test.go
@@ -784,6 +784,10 @@ func TestCalcCellValue(t *testing.T) {
 		"=AVERAGEA(A1)":     "1",
 		"=AVERAGEA(A1:A2)":  "1.5",
 		"=AVERAGEA(D2:F9)":  "12671.375",
+		// BETAINV
+		"=BETAINV(0.2,4,5,0,1)": "0.303225844664082",
+		// BETA.INV
+		"=BETA.INV(0.2,4,5,0,1)": "0.303225844664082",
 		// CHIDIST
 		"=CHIDIST(0.5,3)": "0.918891411654676",
 		"=CHIDIST(8,3)":   "0.0460117056892315",
@@ -859,6 +863,14 @@ func TestCalcCellValue(t *testing.T) {
 		"=FINV(0.5,4,8)":   "0.914645355977072",
 		"=FINV(0.1,79,86)": "1.32646097270444",
 		"=FINV(1,40,5)":    "0",
+		// F.INV.RT
+		"=F.INV.RT(0.2,1,2)":   "3.55555555555555",
+		"=F.INV.RT(0.6,1,2)":   "0.380952380952381",
+		"=F.INV.RT(0.6,2,2)":   "0.666666666666667",
+		"=F.INV.RT(0.6,4,4)":   "0.763454070045235",
+		"=F.INV.RT(0.5,4,8)":   "0.914645355977072",
+		"=F.INV.RT(0.1,79,86)": "1.32646097270444",
+		"=F.INV.RT(1,40,5)":    "0",
 		// NORM.DIST
 		"=NORM.DIST(0.8,1,0.3,TRUE)": "0.252492537546923",
 		"=NORM.DIST(50,40,20,FALSE)": "0.017603266338215",
@@ -2282,6 +2294,32 @@ func TestCalcCellValue(t *testing.T) {
 		"=AVERAGE(H1)": "AVERAGE divide by zero",
 		// AVERAGEA
 		"=AVERAGEA(H1)": "AVERAGEA divide by zero",
+		// BETAINV
+		"=BETAINV()":               "BETAINV requires at least 3 arguments",
+		"=BETAINV(0.2,4,5,0,1,0)":  "BETAINV requires at most 5 arguments",
+		"=BETAINV(\"\",4,5,0,1)":   "strconv.ParseFloat: parsing \"\": invalid syntax",
+		"=BETAINV(0.2,\"\",5,0,1)": "strconv.ParseFloat: parsing \"\": invalid syntax",
+		"=BETAINV(0.2,4,\"\",0,1)": "strconv.ParseFloat: parsing \"\": invalid syntax",
+		"=BETAINV(0.2,4,5,\"\",1)": "strconv.ParseFloat: parsing \"\": invalid syntax",
+		"=BETAINV(0.2,4,5,0,\"\")": "strconv.ParseFloat: parsing \"\": invalid syntax",
+		"=BETAINV(0,4,5,0,1)":      "#NUM!",
+		"=BETAINV(1,4,5,0,1)":      "#NUM!",
+		"=BETAINV(0.2,0,5,0,1)":    "#NUM!",
+		"=BETAINV(0.2,4,0,0,1)":    "#NUM!",
+		"=BETAINV(0.2,4,5,2,2)":    "#NUM!",
+		// BETA.INV
+		"=BETA.INV()":               "BETA.INV requires at least 3 arguments",
+		"=BETA.INV(0.2,4,5,0,1,0)":  "BETA.INV requires at most 5 arguments",
+		"=BETA.INV(\"\",4,5,0,1)":   "strconv.ParseFloat: parsing \"\": invalid syntax",
+		"=BETA.INV(0.2,\"\",5,0,1)": "strconv.ParseFloat: parsing \"\": invalid syntax",
+		"=BETA.INV(0.2,4,\"\",0,1)": "strconv.ParseFloat: parsing \"\": invalid syntax",
+		"=BETA.INV(0.2,4,5,\"\",1)": "strconv.ParseFloat: parsing \"\": invalid syntax",
+		"=BETA.INV(0.2,4,5,0,\"\")": "strconv.ParseFloat: parsing \"\": invalid syntax",
+		"=BETA.INV(0,4,5,0,1)":      "#NUM!",
+		"=BETA.INV(1,4,5,0,1)":      "#NUM!",
+		"=BETA.INV(0.2,0,5,0,1)":    "#NUM!",
+		"=BETA.INV(0.2,4,0,0,1)":    "#NUM!",
+		"=BETA.INV(0.2,4,5,2,2)":    "#NUM!",
 		// AVERAGEIF
 		"=AVERAGEIF()":                      "AVERAGEIF requires at least 2 arguments",
 		"=AVERAGEIF(H1,\"\")":               "AVERAGEIF divide by zero",
@@ -2375,6 +2413,14 @@ func TestCalcCellValue(t *testing.T) {
 		"=FINV(0,1,2)":      "#NUM!",
 		"=FINV(0.2,0.5,2)":  "#NUM!",
 		"=FINV(0.2,1,0.5)":  "#NUM!",
+		// F.INV.RT
+		"=F.INV.RT()":           "F.INV.RT requires 3 arguments",
+		"=F.INV.RT(\"\",1,2)":   "strconv.ParseFloat: parsing \"\": invalid syntax",
+		"=F.INV.RT(0.2,\"\",2)": "strconv.ParseFloat: parsing \"\": invalid syntax",
+		"=F.INV.RT(0.2,1,\"\")": "strconv.ParseFloat: parsing \"\": invalid syntax",
+		"=F.INV.RT(0,1,2)":      "#NUM!",
+		"=F.INV.RT(0.2,0.5,2)":  "#NUM!",
+		"=F.INV.RT(0.2,1,0.5)":  "#NUM!",
 		// NORM.DIST
 		"=NORM.DIST()": "NORM.DIST requires 4 arguments",
 		// NORMDIST
diff --git a/excelize_test.go b/excelize_test.go
index bafd4462477fbec9e20877abb0c57e7f7873ab6a..7d3830495b4862868284e1d75444e87eac6e3520 100644
--- a/excelize_test.go
+++ b/excelize_test.go
@@ -130,14 +130,17 @@ func TestOpenFile(t *testing.T) {
 	// Test boolean write
 	booltest := []struct {
 		value    bool
+		raw      bool
 		expected string
 	}{
-		{false, "0"},
-		{true, "1"},
+		{false, true, "0"},
+		{true, true, "1"},
+		{false, false, "FALSE"},
+		{true, false, "TRUE"},
 	}
 	for _, test := range booltest {
 		assert.NoError(t, f.SetCellValue("Sheet2", "F16", test.value))
-		val, err := f.GetCellValue("Sheet2", "F16")
+		val, err := f.GetCellValue("Sheet2", "F16", Options{RawCellValue: test.raw})
 		assert.NoError(t, err)
 		assert.Equal(t, test.expected, val)
 	}
diff --git a/rows.go b/rows.go
index 0d2490ca651cf8fd14460ce9d8b7b2947a6fbdab..81eaeeb823549db0f02219aca5c2b16b955af389 100644
--- a/rows.go
+++ b/rows.go
@@ -429,6 +429,16 @@ func (c *xlsxC) getValueFrom(f *File, d *xlsxSST, raw bool) (string, error) {
 	f.Lock()
 	defer f.Unlock()
 	switch c.T {
+	case "b":
+		if !raw {
+			if c.V == "1" {
+				return "TRUE", nil
+			}
+			if c.V == "0" {
+				return "FALSE", nil
+			}
+		}
+		return f.formattedValue(c.S, c.V, raw), nil
 	case "s":
 		if c.V != "" {
 			xlsxSI := 0
diff --git a/rows_test.go b/rows_test.go
index 208b2de841377ce6fd584bd8e4f91b33eec32bb9..22b038aa8afb9ebee18d11dfc38bbede54dff0ab 100644
--- a/rows_test.go
+++ b/rows_test.go
@@ -950,6 +950,10 @@ func TestNumberFormats(t *testing.T) {
 	assert.NoError(t, f.Close())
 }
 
+func TestRoundPrecision(t *testing.T) {
+	assert.Equal(t, "text", roundPrecision("text", 0))
+}
+
 func BenchmarkRows(b *testing.B) {
 	f, _ := OpenFile(filepath.Join("test", "Book1.xlsx"))
 	for i := 0; i < b.N; i++ {