Geothmetic Meandian
April 13, 2021
These are straight forward:
(define (agm x y)
(let loop ((a (max x y)) (g (min x y)))
(if (< (/ a g) 1.00000001)
(/ (+ a g) 2)
(loop (/ (+ a g) 2) (sqrt (* a g))))))
(define (gm . xs)
(define (median3 a b c)
(cond ((or (<= b a c) (<= c a b)) a)
((or (<= a b c) (<= c b a)) b)
(else c)))
(let* ((len (length xs))
(a (/ (sum xs) len))
(g (expt (apply * xs) (/ len)))
(m (list-ref (sort < xs) (quotient len 2))))
(let loop ((a a) (g g) (m m))
(if (and (< (/ a g) 1.00000001)
(< (/ a m) 1.00000001)
(< (/ g m) 1.00000001))
(/ (+ a g) 2)
(loop (/ (+ a g m) 3)
(expt (* a g m) 1/3)
(median3 a g m))))))
You can run the programs at http://ideone.com/Z2ZZho.
Fun little drill. I’m not sure about the usefulness of these centrality metrics, but it’s good practice though! Here is my take on this drill using Julia 1.6: https://pastebin.com/sj1X0G0H
(defun arithmetic-geometric-mean (a b) (values (/ (+ a b) 2) (sqrt (* a b)))) (defun converge (epsilon fun a b) (loop :while (< epsilon (abs (- a b))) :do (multiple-value-setq (a b) (funcall fun a b)) :finally (return (values a b)))) (converge 1d-15 (lambda (a b) (format t ";; ~24,21F ~24,21F~%" a b) (arithmetic-geometric-mean a b)) 24d0 6d0) ;; 24.000000000000000000000 6.000000000000000000000 ;; 15.000000000000000000000 12.000000000000000000000 ;; 13.500000000000000000000 13.416407864998739000000 ;; 13.458203932499370000000 13.458139030990985000000 ;; 13.458171481745177000000 13.458171481706055000000 13.458171481725616D0 13.458171481725616D0 (defun geothmetic-meandian (set) (sort (vector (/ (reduce (function +) set :initial-value 0.0d0) (length set)) (expt (reduce (function *) set :initial-value 1.0d0) (/ (length set))) (aref set (truncate (length set) 2))) (function <))) (defun converge-set (epsilon fun set) ;; set is sorted. (loop :while (< epsilon (abs (- (aref set 0) (aref set (1- (length set)))))) :do (setf set (funcall fun set)) :finally (return set))) (converge-set 1d-3 (lambda (set) (format t ";;~{ ~24,21F~}~%" (coerce set 'list)) (geothmetic-meandian set)) #(1 1 2 3 5)) ;; 1.000000000000000000000 1.000000000000000000000 2.000000000000000000000 3.000000000000000000000 5.000000000000000000000 ;; 1.974350485834820000000 2.000000000000000000000 2.400000000000000000000 ;; 2.000000000000000000000 2.116192460544808400000 2.124783495278273000000 ;; 2.079536819479580200000 2.080325318607693800000 2.116192460544808400000 ;; 2.080325318607693800000 2.091948604915222300000 2.092018199544027500000 ;; 2.088090133120960000000 2.088097374355647700000 2.091948604915222300000 ;; 2.088097374355647700000 2.089377914218486500000 2.089378704130609800000 #(2.088951243615992D0 2.0889513309015815D0 2.0893779142184865D0)You can run it at: https://ideone.com/8s3uh4
Here’s a solution in Python, utilizing the
statisticsmodule added to the standard library in Python 3.4.from statistics import geometric_mean, mean, median EPSILON = 1e-10 def geometric_meandian(num1, num2): while abs(num1 - num2) > EPSILON: num1, num2 = mean((num1, num2)), geometric_mean((num1, num2)) return num1 def geothmetic_meandian(*nums): while max(nums) - min(nums) > EPSILON: nums = mean(nums), geometric_mean(nums), median(nums) return nums[0] print(geometric_meandian(24, 6)) print(geothmetic_meandian(1, 1, 2, 3, 5))Output:
Here’s another Python solution, which doesn’t rely on the
statisticsmodule.import math EPSILON = 1e-10 def geometric_meandian(num1, num2): while abs(num1 - num2) > EPSILON: num1, num2 = num1 / 2 + num2 / 2, math.sqrt(num1 * num2) return num1 def geothmetic_meandian(*nums): nums = sorted(nums) while nums[-1] - nums[0] > EPSILON: amean = sum(x / len(nums) for x in nums) gmean = math.prod(nums) ** (1 / len(nums)) median = nums[len(nums) // 2] if len(nums) & 1 == 0: median = median / 2 + nums[len(nums) // 2 - 1] / 2 nums = sorted([amean, gmean, median]) return nums[0] print(geometric_meandian(24, 6)) print(geothmetic_meandian(1, 1, 2, 3, 5))Output:
I have written the code of the following problem in Fortran using subroutine, so that it can be used for finding the solution of other inputted numbers also
<a href=”https://pastebin.com/ZQ7inPa4”>
Here’s a Haskell version. I used this as an excuse to play around with the foldl library.
import Control.Applicative (liftA2) import qualified Control.Foldl as F import Data.Bool (bool) import Data.List (sort) import Data.List.NonEmpty (NonEmpty, fromList, toList) import Data.Maybe (fromJust) {- ----------------------- Arithmetic-Geometric Mean ------------------------ -} agm :: (Floating a, Ord a) => a -> (a, a) -> a agm tol = fst . while ((> tol) . diff) step where step (a, g) = ((a + g) / 2, sqrt (a * g)) diff (x, y) = abs (x - y) {- --------------------------------- Median --------------------------------- -} medianF :: (Fractional a, Ord a) => F.Fold a a medianF = F.fold F.mean . middle <$> sortF where sortF = (,) <$> F.length <*> (sort <$> F.list) middle (n, xs) = take (bool 1 2 (even n)) $ drop ((n - 1) `div` 2) xs {- ----------------------------- Geometric Mean ----------------------------- -} gmeanF :: Floating a => F.Fold a a gmeanF = nthRoot <$> F.product <*> F.length where nthRoot x n = x ** (1 / fromIntegral n) {- -------------------------- Geothmetic Meandian --------------------------- -} gmm :: (Floating a, Ord a) => a -> NonEmpty a -> a gmm tol = head . gmms . toList where gmms = while ((> tol) . maxDiff) (F.fold gmmF) gmmF = sequenceA [F.mean, gmeanF, medianF] maxDiff = fromJust . F.fold (liftA2 (-) <$> F.maximum <*> F.minimum) {- -------------------------------------------------------------------------- -} while :: (a -> Bool) -> (a -> a) -> a -> a while p f = head . dropWhile p . drop 1 . iterate f test :: Show a => String -> (a -> Double) -> a -> IO () test name f arg = putStrLn $ name ++ ": " ++ show arg ++ " => " ++ show (f arg) main :: IO () main = do let tol = 0.00001 xy = (24, 6) ys = fromList [1, 1, 2, 3, 5] test "AGM" (agm tol) xy test "GMM" (gmm tol) ys